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ABSTRACT 


Power algorithms for several common statistical tests are presented. In 
experimental design of operational tests and evaluations the selection of design 
parameters so as to attain an experiment with desired power is a difficult and 
important problem. 

An interactive computer program is presented which uses the power algorithms 
for several tests and creates graphical presentations which can be used to assist 
decision makers in statistical design. Several common tests and associated parameters 


(such as sample size. types and levels of treatments, and alpha-level) are examined. 
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I. INTRODUCTION 


A. DESCRIPTION OF THE PROBLEM 

A statistical design is a plan according to which an experiment is patterned. It 
provides the basis upon which appropriate statistical tests and inferences can be made 
after the experiment has been performed. The selection of the experimental design to 
be used in a given situation is extremely important because it plays a predominant role 
in the efficiency of the experiment, the precision with which the objectives are met, and 
the total effort (and cost) expended upon the experiment. 

The concern of this thesis is the development of power algorithms for several 
common tests, including the F-test and t-test. The experimenter wants a test to 
achieve the correct decision with as high probability as possible in practical operational 


tests and evaluations given certain conditions hold. 


B. SCOPE OF THE THESIS 

An interactive computer program is developed which uses power algonthms for 
several common tests. The program can be used by decision makers who may not 
have deep knowledge of statistics or may be unfamiliar with the details of the design of 
experiments. The output helps the decision maker design experiments appropriate for 
testing hypotheses with common statistical tests. The experiment designer can use the 
program to evaluate the effect upon power of varying underlying design parameters 
such as sample size. 

Because the underlying statistical noncentral distributions used in computing 
power do not have tables, approximation methods are used for computing the 
noncentral CDFs. Approximation methods were developed which are efficient, 


requiring reduced computer CPU time. 


C. BACKGROUND 
A hypothesis is a statement about the values of the parameters of a probability 
distribution. For example, suppose we think that the mean yield of a chemical process 
is more than 94.5 percent. Hypotheses to test this statement might be expressed 
formally as 
Ho: S 94.5 
Ha:p > 94.5 


The statement Ho: w S 94.5 1s called the null hypothesis, and Ha: pt > 94.5 is called 
the alternative hypothesis. The value of the mean specified in the null hypothesis 
might typically be determined in one of two ways. It may be the result of some theory 
or model regarding the process under study, or it may be the result of contractual 
specifications. 

To test a hypothesis one usually devises a procedure for taking a random sample, 
computes an appropriate test statistic, and then rejects or declines to reject the null 
hypothesis Ho, depending on the outcome on the test statistic. Part of this procedure 
is the specification of the set of values for the test statistic which lead to rejection of 
Ho. This set of values is called the critical region or rejection region for the test. 

Two kinds of errors may be committed when testing hypotheses. If the null 
hypothesis is rejected when it is true, then a type I error has occurred. If the null 
hypothesis is not rejected when it is false, then a type II error has been made. The 
probabilities of these two errors are given special symbols: 

@ =P(type I error)= P{reject Ho | Ho 1s true) 
{f} = P(type II error)= P(fail to reject Ho | Ho is false) 
Generally, the determination of B requires the use of a noncentral sampling 


2 , and noncentral t) for which tables are 


distribution (e.g., noncentral F, noncentral ¥ 
not readily available. The noncentral sampling distribution depends on a parameter 
called the noncentrality parameter. The noncentrality parameter is usually a measure 
of the distance {in some sense) between the values of the parameter under the null and 
alternate hypotheses. Thus, type I] error rates depend upon parameters of the 
distribution of the test statistic under the alternative hypothesis. The power of the test 
for a spevified value of the hypothesized parameter is the probability the test would 
reject Ho when Ho is false, or 1 — fi. 

In addition, whatever the test procedure is, rejecting Ho when Ho is not true is 
usually something we want a test to achieve with as high a probability as possible. 
Therefore we want the power of the test, for a given value of non-centrality parameter, 
to be high. Power algorithms are used to assist in the “best” selection of parameters 
for a statistical design. There do not exist simple closed form equations for computing 
the required probabilities from the common noncentral sampling distributions. 
Consequently, one must use numerical approximations to determine power. In this 
thesis we use various approximation methods to compute the cumulative probabilities 


2 


for the noncentral F, noncentral y “, and noncentral t distributions. 
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Il. NONCENTRAL CHI-SQUARE,T,F DISTRIBUTIONS 


A. NONCENTRAL CHI-SQUARE DISTRIBUTION 

If W1, W2, .... Wn are independently distributed as NC) ele 2). 1, then 
yy Wi2 has a distribution known as a noncentral 4? (n, % ) where n is the number of 
degrees of freedom and’ = Y pt 2 is the noncentrality parameter. When! ; = Ht 5 
ae = ft, = 0, then dh = QO, and the noncentral 42 (n,O) reduces to the usual central 


42 (n) with n degrees of freedom. The cumulative distribution function of 4? (n,Xr ) is, 


F(%;) 9,4) = 


| 
=v 
TaaN 
% 
cos 
s 
>” 
in 
~ 
—! 
v 








while F(x:n, 4) = 0 for x S$ 0 


It is possible to express F(x;n,A4) for x>0, in an easily remembered form as a 


y 


weighted sum of central y~ CDF’s with weights equal to the probabilities of a Poisson 


distribution with expected value 4/2. That is, 





oo Ly) ty 2 
Berths © (8° | PU Xeveap 4) (2.2) 


7! 
Thus a 47 (n, 4) variable can be regarded as a mixture of central 42 variables. This 
interpretation is often useful in deriving the distribution of functions of noncentral 42 
random variables. 

The probability density function can, similarly, be expressed as a mixture of 


central 4? probability density functions 


= Nath ge? 2 
fia) = OE, where (= 


ea). S. a (2.3) 
j= [Ty rt) 27 eta 2= 


el 


The mean and variance of the distribution are 
E(X)=n+2A and Var(X)=2n+ 8A [Ref. I:p. 130]. 
NONCENTRAL F DISTRIBUTION 
If Y,; and Y> are independent and Y, 1s x | 2 (n, ,A ) and Y> 18 X 2 (n, ) then 


Va (ea OG ae, (2.4) 


is distributed as F’(n, nN, A), the noncentral F distribution with ny and n, degrees of 
freedom and noncentrality parameter i. Its density function is given by 





co Non pt (a+ Be +k] yk-t 
pW) = Ce * VzZoO 25 
=e h (2: +K) (notny)® , ( " ) 
where a n> fo pe 


c= 
(BS) (natn) tae 


when A =0 this reduces to the density function of the central F(n,,n, ) distribution. 
The mean and variance of the distribution are 


1+ 
N.-2 Ny 





E(v)= —4 aa 


fOr n= 2 


and 


2 a 
Var(V) = iu (M+ 2A) a ten (for ho = 4) 
Ny (M2) CNar2)(Na- 4) Na- |} 


When A =0, these reduce to the mean and variance of the central F(n,,n,) 
distribution. Derivation of (2.5) is shown in [Ref. I:p. 189]. 





The cumulative distribution of V can be expressed in terms of an infinite series of 
multiples of incomplete beta functions, as follows: 


= (4A? =), 
REV sfJ=E (a o)- Tn (Ane, 2) 


Natnit = Ce 


where I wy ( nj, me is the incomplete beta function. 
= Wor 


In this thesis, we use Paulson’s Approximation from Severo and Zelen [Ref. 2]. 


Their approximation 1s as follows: 


PoV Sf,) ~ P(x) (2.7) 


w]e 


where a 
oe ae ae a = a ile | 


(= 


, and 


® is the standard Normal CDF. 

The following table shows some values obtained using this approximation, 
together with exact values of Pr(V = fo ). It can be seen in TABLE 1 that this 
approximation gives about 2 decimal place accuracy. Thus, this approximation 
provides adequate accuracy for computing the power of F-tests [Ref. 3:p. 202]. 

The program we used in the computation of noncentral F CDF values 1s shown 


in Appendix D. 


C. NONCENTRAL T DISTRIBUTION 


The ratio 
seyleigre 20)) x a(n) (2.8) 


where U and Y are independent random variables distributed as standard normal 
(N(0,1)) and x 2 with n degrees of freedom ,respectively, is said to have a noncentral t 
distribution with n degrees of freedom and noncentrality parameter A. Sometimes 2 
(orc ven 5h), rather than A, is termed the noncentrality parameter. If 4 is equal to 
zero, the distribution is central t with n degrees of freedom. 

The cumulative distribution of T 1s [Ref. I:p. 201] 
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PrT St) = 


Tele al 
CDF OF A NONCENTRAL F DISTRIBUTION 


CXaegl 


approx 


fo 


n2 


nl 





20 





Abamowitz and 


ggested by 


we use an approximation method su 


Stegun [Ref. 4:p. 949]. 


3 


In this thesis 


10) 


2 


( 


t) ~ @ (x) 


= 


Pr 


and @ is the standard Normal CDF. 


where x 
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The following table shows some values obtained using this approximation, 
Pee@ciner with ckxact values Of Prt = to ). It can be seen in Table 2 that this 
approximation gives about 2 decimal place accuracy. Thus, this approximation 
provides adequate accuracy for computing the power of t-tests. 


The program we used in the computation of noncentral t CDF values is shown 
Appendix E. 


TAB RE 2 
CDF OF A NONCENTRAL T DISTRIBUTION 





n Xd to approx exact 
4 0.112 4 0.017 0.01 
7.944 4 O95 0.99 
2.24 10 0.026 0.01 
18.494 10 eS 99 
5.78 79 0.012 0.01 
one 79 Oeil 0.99 
16.1 4 0.014 O] 
Doe 0.992 ae 


Il. POWER OF THE F-TEST 


A. INTRODUCTION 

Under appropriate conditions, the best test for testing equality of several means 
is the analysis of variance test. Analysis of variance has a wide application. It is one 
of the most useful techniques in the field of statistical inference. As in any 
hypothesis-testing situation, the power of the F test is of interest to the experimenter. 

In this chapter we will discuss the power of F tests and provide an example. To 
give an overall evaluation of the power of F tests in the analysis of vanance, we may 
use power curves. An important use of the power curve is to guide the experimenter in 
selecting the sample size (number of replicates) so that the design will be sufficiently 
sensitive to important potential differences in the treatments. We will consider power 


curves for one-way and multi-way analyses of variance (ANOVA’S). 


B. THE ONE-WAY CLASSIFICATION ANALYSIS OF VARIANCE 

Suppose we have k different levels of a single factor that we wish to compare. 
The different levels of the factor are often called treatments. The observed responses 
from each of k treatments is a random sample on a random variable. The data would 


appear as in Figure 3.1. 


Observation 


Treatment 





Figure 3.1 Typical Data for One-way Classification ANOVA. 


We will find it useful to describe the observations by the linear statistical model: 


Be i eG (3.1) 


yth observation, {f is a parameter common to all treatments called 


where Yij is the (yj 
h h 


treatment called the i” treatment 


2), 


the overall mean, Tt; is a parameter peculiar to the i! 


effect, and ei is a random error component, assumed to be IID N (0, 6 


The indices used are: 


eee) — tie treatmients, 1 = 1,2....,k 


ej Hie me plication per treatment, j — 1,2,...,n 

The objective in the ANOVA is to test appropriate hypotheses about the 
treatment effects. The variance o* is assumed constant for all levels of the factor. 
This model is called the one-way classification analysis of variance because only one 
factor 1s investigated. 

We are interested in testing the equality of the k treatment effects, so the 


appropriate hypotheses are 


Ho: ty = Ue) er = ="C 


[ae ¢; for some 1,} . 


That is, if the nuil hypothesis is true, then each observation 1s made up of the mean 
+ tT plus a realization of the random error Eji. 
The results of the ANOVA procedure is summarized in Table 3. 
We may use the expected values of the mean squares to verify that Fo (in Table 
3) is an appropriate test statistic for Ho. From the expected mean square we see that, 


in general, MSE is an unbiased estimator of G2 


. Also, under the null hypothesis, MSt 
is an unbiased estimator of 6%. However, if the null hypothesis is false, then the 
expected value of MSt 1s G7 +(n ¥ (t;)7)/(k-1). diese pected valde of the mean 
square error is 67. | 

Therefore, under the alternate hypothesis the expected value of the numerator of 
the test statistic (Fo) is greater than the expected value of the denominator and we 
would reject Ho with values of the test statistic which are too large. That is, we would 


reject Hopi: 


By Fa(k-1)(N-k) 
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TABLE 3 
THE ANALYSIS OF VARIANCE FOR A ONE-WAY CLASSIFICATION 


Source of Sum of Degrees of Mean 
Variation Squares Freedom Square 


Between Sse k=] MSt MSt/MSE 
Treatments 


Error(within Soe MSE 
treatments ) 


Total Soi 





where the numerator and denominator df are k-| and k(n-1)=N-k, and @ 1s the type I 
Crone 


The power of the test 1s: 
1-8 = rer > F a(k-1)(N-k) | Ho is false) (3.2) 


To evaluate the B in Equation 3.2 we need to know the distribution of the test statistic 
re if the null hypothesis is false. It can be shown that, if Ho 1s false, the statistic Ea 
has the noncentral F distribution with k-I and k(n-1) degrees of freedom and 


noncentrality parameter i, given by 


nY(t;) 





62 


The noncentrality can be interpreted as the squared standardized distance 
between the origin and (Tj,T5, ... ,t}). The ratio y (tT 2 | o? is called the squared 


2 


standardized distance. If only an estimate of o% is available, one may replace 6“ with 


the estimate [Ref. 5:p. 34]. 


C. THE MULTI-WAY CLASSIFICATION ANALYSIS OF VARIANCE 
Many experiments require a study of the effects of two or more factors. It can 
oe shown that, under certain conditions, factorial arrangements are the most eflicient 


designs for this type of analysis. 
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One of the simplest factorial experiments involves only two factors or sets of 
treatments say factor A and factor B. Suppose there are a levels of factors A and b 
levels of factor B, and these are arranged in a 2-way factorial design; that is, each 
replication of the experiment contains all ab treatment combinations. Assume there 
are n replications of the experiment, and let Yajk represent the observation taken under 
the i'® Jevel of factor A and the jth level of factor B in the kth replication. 

The data can be summarized as shown in Figure 3.2. The order in which the abn 


observations are taken 1s selected at random. 


Factor B 
| | 1 | 2 | | b | 
ee | ere | | 
Pee ime 24>')22 | | Vapi Yip2 | 
ME \ =. ee eT 
| ee || ee | | 
Pe ON eel 201 3N200) || | Yop1>Yap2 | 
me ep ee eee Youn | 
morc eee | | 
A ee. | | | | 
Lgeciines | | | | | 
Peer | | | | | 
| Le |e ee | 
ee io aaa! 'ag]>'aae | | Yab1>%ab2 | 
| a | ’ ’ On | ) ) ean | | ’ ) lean | 
| | | | 


Figure 3.2 Typical Data notation for a Two-way Classification.. 


The observations may be described by the linear model: 
comet tt Pp te, ea) 
where 


i 


e jt = the overall mean effect; 


¢ t; = the true effect of 1 th jevel of factor A,1 Pee 
e f ‘i. the true effect of j CD evel of factor Bie wl 2, ce 0: 
e th a the effect of the interaction between t ; and 3 and 
° fi, = a random error component, assumed to be IID N (0,6 v yy 
Both factors are assumed to be fixed. It is usually assumed that the treatment effects 
are defined as deviations from the overall mean, so ) t ; = 0 and » fp io 0. 
Similarly, the interaction effects are fixed and usally defined so that y (tp) jn oy 
Where the summation 1s over either 1 or j. Since there are n replicates of the 
experiment, there are a total of abn observations. 

We are interested in testing various hypothesis about the parameters in equation 
3.3. An appropriate hypothesis testing procedure would again be analysis of variance. 
More specifically, as we are considering two controllable sources of variation (A and 
B), the procedure is called the two-way classification analysis variance. 

In Order to test the hypothesis om... ae for 1=1,2,...,a (nO TOW factor 


effects), Ho; B a QO for }=1,2,...,6 (no column factor effects), and Ho ;( tf ) ij = 0 


(no interaction effects), we can express the total sum of square as: 
SST = SSA“ SSB > SSAB 4 oe (3.4) 


Here, SSA is a sum of squares due to “rows” or factor A, SSB is a sum of squares due 
to “columns” or factor B, SSAB is a sum of squares due to the interaction between A 
and B, and SSE is a sum of squares due to error. The degrees of freedom associated 
with each sum of squares are shown in Figure 3.3. 

If we assume € ijk are ITD N(0,67 ) and apply Cochran's theorem (Theorem 3-1) 
under the null hypothesis of no effects, each sum of squares on the right-hand side of 
Table 4 when divided by o* is distributed as x 2 with the indicated number of degrees 


of freedom, and these statistics are independent. 


Theorem 3-1 (COCHRAN). Let Z ; be IID N(0,1) for 1=1,2,....v and suppose 
WA 2 = Q; + Q + Q; +...+Q, where s S v, and the quadratic form Q; has 
v, degrees of freedomi(1 [2.25 ea ten Q), Q>5 .., Q, are independent 
chi-square random variable with Vv; V>.., Vv. degrees of freedom, respectively, if 
and only if 

Se ae “7 a 5 Se rare eaie Vo. 
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Effect Degrees of freedom 

a 
B Sai 

AB interaction (aaa 1) 


Enron ab(n-1) 





Total abn-1 | 


Figure 3.3 Table of degrees of freedom with sums of squares.. 


Assuming that factors A and B are fixed, the expected values of the 


mean squares are: 


E(MSA) = o% +( bn Yt; 7)/(a-1); (3.5) 
E(MSB) = o” +( bn YB; 7 )/(b-1) ; (3.6) 
E(MSAB) = 6% + (n Oa) ' 2 y/ (a-1)(b-1) ; and (3.7) 
E(MSE) = o%. (3.8) 


Therefore, to test the hypothesis Ho: T eae ae i= 1,2,...,a (no row factor effects), and 
Ho ; B 


interaction effects), we would divide the corresponding mean square by the mean 


= 0; j=1,2,...,5 ( no column factor effects), and Ho: ( t B ) a 0 (no 


square error. Under the null hypothesis of no effect, this ratio will follow an F 
distribution with appropriate numerator degrees of freedom and ab(n-1) denominator 
degrees of freedom, and the critical region will be located in the upper tail. The test 
procedure is usually summarized in an analysis of variance table, such as shown in 
Table 4. 


Z| 


TABLE 4 
THE ANALYSIS OF VARIANCE FOR A TWO-WAY CLASSIFICATION 





Source of Sum of Degrees of Mean 

Variation Squares Freedom Square 

A treatments SSA a-1 MSA MSA/MSE 
B treatments SSB Pe MSB MSB/MSE 
Interaction SSAB (a-1))( 5-1) MSAB MSAB/MSE 
Erron oe ab(n-1) MSE 
Tota] Soi abn-1 


To the compute the power of tests in two-way ANOVAs. the procedure is the 


same as in the one-way case. Relationships among A, the numerator degrees of 


freedom and the denominator degrees of freedom are shown in Table 5. 


In a similar way one can expand to the general multi-way ANOVA procedure 


[Ref. 5:p. 124]. 


D. 


THE ALGORITHMS AND FLOWCHART 


A program for computing a power table and power curve in general ANOVA’s 


(fixed model) is shown in Appendix F. The power of one-way ANOVA’s are solved 


with these programs, using the following the sequence: 


I) 
2) 
3) 
4) 


5) 


Determine what variable to include, such as sample size vs power, number of 
treatments vs power, @-level vs power, or noncentrality vs power. 


Detgmaamne the values of the maximum, minimum, and increment of the variable 
chosen. 


Given input data, compute the critical value, using the inverse of the central 
F-distribution (Appendix A). 


nua the power value in each case; that is, the CDF of the noncentral 
F-distribution iene: ips 


Print the power table or the power curve. 


ee 


TABLE 5 
NONCENTRALITY PARAMETERS FOR POWER IN A TWO-WAY ANOVA 


Numerator Denominator 
Factor iz DOF DOF 
2 
Da t. 
A pnd = a-l ab(n-1) 
Oo 
Z 
an 
B ae Dal ab(n-1) 
ny tBs 5° 
EiesS 3 — = Goo) ab(n-1) 
O 


Multi-way ANOVAs have the same algorithms but multi-way ANOVA may 
include tests of interaction effects. The program considers only up to three-way 
interaction effects. We can explain how to compute degrees of freedom of error term 
in m-way ANOVA’s, assuming the balanced case, as follows: Total degrees of freedom 
is DOF(Total)=n » k(i) -1, where k(1) 1s the number of levels of the MD factor Pixel) |: 

(1) If the model has only main effects without any interaction effects, 
DOF 1(Error)= DOF(Total) - ¥° (k(i)-1). 

(2) If the model has several factors and only 2-way interaction effects, 
DOF2(Error)= DOF (Error) - ¥) ¥ (k(i)-1)(k(j)-1). 

(3) IF the model has several factors and only 3-way interaction effects, 
DOF3(Error)= DOF 1(Error) - © YY (k(i)-1)(KG)- 1)(K(K)- 1). 

(4) If the model has several factors and 2-way and 3-way interaction 
effects. 
DOF4(Error) = DOF3(Error) - YY (k(i)-1)(k(j)-1). 

If the model has more than 3-way interaction effects, then the user must modify 
the degrees of freedom of the error terms accordingly. A flowchart is shown in Figure 
3.4. 


i 


START 


Determine one-way or / multi-way 
multi-way ANOV ? 
one-way How many ways (M)? 


For i-factor level, 
How many level (MF)? 












Determine what to include? 


— nee ee = 
—_— =~ 






sample size vS power Determine how to 


. alpha-level vs power many way interactinn? 


. Moncentrality vs power 








Compute corresponding to 


degrees of freedom 





Determine what to include? 


Compute the critical value 
and power 


a 


# of replicates vS power 
noncentrality vs power 


alpha-level vS power 






Plot power curve . 
Compute the critical 


value and power 


Analyze again? 


STOP 


Figure 3.4 System flowchart for F-test power. 
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E. EXAMPLES OF F-TESTS 

1. Example of One-way ANOVA test (number of replicates vs power) 

a. Scenario 
A manufacturer suspects that the batches of raw material furnished by his 
supplier differ significantly in calcium content. There are a large number of batches 
currently in the warehouse. Five of those are randomly selected for study. A chemist 
wants to Know the appropriate sample size per each batch, in order to test 
OMG = U5 Gon Ty 

against 


ent T; for some 1,] 


The chemist would like to know how many replicates to run if it 1s important to reject 
Ho with probability at least 0.9 when the standardized distance square ( ) t 2 i) 2 ) 
is 2 and @ = .05. Thus he would like to Know what the power of the F-test is for a 
range of possible replicates. He decides to check replicates from 2 to 12 in increments 
oF 1. 
b. Inputs 
1) Select one-way ANOVA test (the number of replicates vs power). 
2) Ihe number of treatment = 5. 
3) @-level = .05 
4) Standardized distance square = 2. 
5) M>->ximum replicates is 12, minimum replicates is 2, increment 1s 1. 
c. Start program 
Do you want to analyze one way ANOVA(y/n)? 
y 
Do you want to plot n (# of observation per treatment vs power)(y/n)? 


iy 
Number of k( # of treatment)? 
t 
5 
Alpha- level? 
? 
OS 


Standardized squared distance value ? 
; 
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2 
Maximum n value (the maximum value on X-axis)? 
? 
12 
Minimum n value (the minimum value on X-axis)? 
(condition :N must be more than 2) 

? 
2 
Increment n value ? 
? 
l 

d. Output 

The output is presented in Table 6 


TABLE 6 
OUTPUT OF THE ONE-WAY ANOVA EXAMPLE (REPLICATES VS POWER) 


a 


DOP IT MEDOr2 ae F-INVERSE POWER 


a) 
bo 
\O 
S 
OO 
=) 


SS 
WBOAnNON D~IOOO 
BOaAN SWODA~TIA 


Om~=—WAW LP — On) 
NW Da NGI Qn 


a 
Sa So oe 
QQ Qoe Oo 


4 
4 
4 
4 
d 
4 
4 
4 
4 
4 
4 


Cacn BS BWUNRND—— Wn 
CAC Ui Cn Cn ic 
PEEP EPP EELS 
OOD Ss 2 a SS Sa Sas SG 
CAC CA Cain n Cnn nnn 
NNNNPNNNWer 
CARINACAOQN~)1QO\O = 
BNE OHDBNOODL 
SOOSOS99999 
O'WD'\O'\0 00 ~~} ~] On BG 
COD LP —=AODo hen 
=~ JONG BU Gnrodw 
mt JONI DO Ba CNGI~J WK ~] 
COD COW DI NI~)J On BW 





# Of replicate POWER 


NN= 2. POWER= _ 0.15373 
NN= __ 3. POWER= _ 0.30024 
NN= 4. POWER= _ 0.44935 
NN=_ 5. POWER = "0 scon 7 
NN= 6. POWER=_ 0.70332 


26 


HF I5118 
POWER 


0.56746 


Oelo373 


NN= = 7. POWER= 0.79562 
INN= 8. POWER=_ 0.86449 
NN=_ 9. POWER= _ 0.91328 
NN= _ 10. POWER= _ 0.94631 
NN= 11. POWER= _ 0.96776 
NN= 12. POWER= _ 0.98118 
ta-fercccee-- foeee--<--- pesecer---- -ecen----- ferere---- foeoX*----- % am 
| i | 
| * | 
| | 
| - | 
| | 
| | 
| | 
| | 
| | 
+ * + 
| | 
| | 
| : | 
| | 
| | 
| . | 
| | 
| | 
| | 
$_X* ence enn--- +ecece---- +c ece---- +eecece---- teecer---- teeece---- +—+ 
2.000 7.000 # OF REPLICATES 12.00 
X-SCALE: *-"= 0.125E+00 UNITS 
Y-SCALE: *|"= 0.138E-01 UNITS 


a 


2. Example of one-way ANOVA test (noncentrality vs power) 
a. Scenario 
Five brands of batteries are under study. It is suspected that the life (in 

weeks) of the five brands is different. Five batteries of each brand are tested. A 
manufacturer wants to Know the power as a function of the standardized squared 
distances. That 1s, would like to Know what the power of the F-test is for a range of 
possible standardized distances square. He decides to check standardized squared 
distances from 1 to 5in .2 increment, where n=5 and @ = .05. 

b. Inputs 

1) Select one-way ANOVA test (the standardized distance square vs power). 


2) Number of treatment = 5. 
3) Number of replicates = 5. 
4) a-level = .05. 


5) Maximum standardized distance square is 5 and minimum standardized distance 
square is | and the increment ts .2. 


c. Start program 
Start program: 
Do you want to analyze one-way ANOVA(y/n)? 
y 
Do you want to plot n(# of observation per treatment) vs power(y/n)? 
n 
Do you want to plot alpha-level vs power(y/n)? 
n 
Do you want to plot noncentrality vs power(y/n)? 
y 
# of observations per treatment (n= ) ? 
? 
» 


# of treatment (k= )? 
? 
5 
Alpha-level ? 
? 
OS 


Maximum standardized squared distance value range ? 
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g 


5 


Minimum standardized squared distance value range ? 


? 


1. 


Increment standardized squared distance value range ? 


: 


d. Output 


The screen output ( Table 7 ) is as follows: 


O05) Boe 


INCENTRALITY VS 


NO 


( 


NOVA EXAMPLE 


POWER) 


WAY A 


OUTPUT OF THE ONE 


DOF2 a 


DOF 


CONST NO SF DIMM NAN OCOP DD) St ENOVNO 


Oo2OO OOOO 002002 
OOOO O00 OO O92 
NOAAMNTOOr-ONnG— 


NoNoNoNeNoNoNoNeNoNoNoNoNoNoNoNoNoNoNoNoNe 
PEE EEE eee eee 
WWD O OOD DODO OOD OOO OO OOO 


ee Sf Sl ee ee ee ee eed eee Sl ee eee ee ee ee et 


ANDANARAAHAARLAANAHAONAHAONONAON 
CQAAANANANANAANANA’A NA ANNA|ANNAN 


DAWN WDD WY YW YY) 
BOOOO00000004000C2 004000 


Tease sg ggg ese ITI TIT TST 
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NONCENTRALITY 


NN= 
NN= 
NN= 
SS 
NN= 
NN= 
iN 
NN = 
NN = 
NN= 
NN= 
NN= 
NN= 
NN= 
NN= 
NN = 
NN = 
NN= 


5.0 

6.0 

7.0 

8.0 

9.0 
10.0 
11.0 
1230 
13.0 
14.0 
LO.) 
16.0 
17.0 
18.0 
oy 
20.0 
2a 
22.0 
23.0 
24.0 
Zo 


POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
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POWER 


0.30382 
0.36314 
0.42203 
0.47946 
0.53461 
0.58677 
0.63550 
0.68050 
UR os 
0.75885 
O72 25 
0.82192 
0.84812 
0.87108 
0.89107 
0.90835 
O77 2528 
Oe. 
0.94672 
0.95586 
0.96356 


0. 96356 
POWER 


0. 63369 


oe s0382 


$<=4+==---=--- $oecemceee-- fescece---- $rernne--- $oecce---- +ece--- % ae oe Kamp 
| ca. | 
| a | 
| | 
| | 
| . | 
| * | 
| . | 
| | 
| : | 
~ * + 
| ‘ | 
| | 
| | 
| | 
| : | 
| . | 
| | 
| | 
| | 
$—_* eee === +-2=------- +22e------ +--------- $oece----- +een------ +=—+ 
5. 000 15.00 NONCENTRALITY 25.00 
X-SCALE = 0.250E+00 UNITS 
Y-SCALE = 0.110E-01 UNITS 
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IV. POWER OF THE T-TEST 


A. INTRODUCTION 

[In testing the hypothesis that the mean of a normal distribution is equal to Hy 
when the standard deviation is unknown, the t statistic may be used. This test-statistic 
is a function of the sample mean, the sample standard deviation, and the sample size. 

The t-test has a very wide application in research areas. For example, we may 
wish to know whether product A 1s better than product B; or whether the outputs of 
machines C and D form a homogeneous mass of product; and so forth. To answer 
such questions we can employ the t-test which is one of the most useful techniques in 
the field of statistical inference. 

In such hypothesis-testing situations the power of t-tests is usually of interest to 
the experimenter. We may use the power curve to evaluate the power of the t-test. 
The required sample size might be determined by referring to a set of power curves. 


We will consider the power curves for the one-sample and two-sample cases. [Ref. 10] 


B. ONE-SAMPLE T-TEST 

We assume that a random sample, XI, X2, X3, .... Xn of size n is taken from a 
normal population with mean pw and standard deviation 6. Suppose we wish to test the 
hypothesis that f# S fp against the alternative hypothesis that f > ply . The 
procedure is to reject the hypothesis Ho: S py foc Hy ) yo) See q » andto 
accept Ho if otherwise where X and S are the sample mean and sample standard 


th 


deviation and where t g is the 1 — @”” percentage point of the Student t-distribution 


with n-1 degrees of freedom, i.e., 
P( Studentt > t,) = @ (4.1) 
With critical region ( ty » © ), the hypothesis f = fp is rejected with probability 


a when Ho is true. If the mean of the normal distribution is Ht, > Hy , then the 


power of the above test is 


P( noncentralt > tg |6 = ((H,- My) Vn)/ 6) = 1-8 (42) 


oe 


where 6 is the noncentrality parameter for the noncentral t-distribution with n-1 
degrees of freedom [Ref. 6:p .320]. 


C. TWO-SAMPLE T-TEST 

Suppose we are given two random samples from two normal populations X1, X2, 
meemeanad yi, 12, ..., ym, where the As and Y § are independent with mean values 
My = By and common unknown variances o, 2 = oy 2 = gt Suppose we wish to 
test the null hypothesis Ho: py S Hy against the alternative hypothesis Ha that py 
> Hy. The procedure is to reject the hypothesis Ho if (Ox vs feet Sp - (Gi) acl ma) 
> t 


where t ,, iS a critical value of the Student t-distribution based on m+tn- 2 


a a 
degrees of freedom. In this case take 6 = ( H — Hy, y's 7 aGivn + 1/m) so the 


probability of rejecting the hypothesis when ph, > Hy, is equal to 
P( noncentralt > ty, |6) = I-B (4.3) 


where 6 is the noncentrality parameter for the noncentral t-distribution with m+n-2 
Meerees of freedom (Ref. 6:p .321]. 


D. YWO-SIDED T-TEST 
For tests where the alternative hypothesis specifies that MH # [ly in the test of 


section B or H, * ft, in the test of section C, where the direction can be either up or 


yi 
down, Le, H > Hy or ft < Hy and fl, > Hy or Wy < Hy, we have what is known as a 
two-sided test. In this case the rejection procedure for the test of section B specifies 
that Ho is rejected if either ((X- py )/ n)/ S > t g/> or (X- wy) Vn)/S < -tg/s 


and the power of the test 1s 
P( noncentral t > t g/>|6) + P( noncentral t < -t,/>{0)= 1-8 (4.4) 


where 6 = (fi, -My) V n/ [Ref. 6:p .323}. 


E. THE ALGORITHMS AND FLOWCHART 

Listings of programs for generating the power table and curve for t-tests are 
given in the Appendix E. 

The following steps are used to compute the power: 


1) Determining whether one-sample or two-sample t-tests should be used. 
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3) 
4) 


5) 


6) 


After determining the test, give the input data. (n,@, Hy), © ) 
Determine whether the one-sided or two-sided test should be used. 


Compute the critical region (inverse t-values). (The computation method is 
provided in Appendix B 


Compute, the power corresponding to the chosen user option. (The 
computation method 1s discussed in Chapter IT) 


Plot the power corresponding to each above option. 


A flowchart is shown in Figure 4.la 


EXAMPLES OF T-TESTS 


1. Example of one sample one-sided t-test 


a. Scenario 


The supervisor in an electronic instrument plant is concerned with the 


repair time of an electronic equipment. Standards specify that the mean repair time 


must be at most 40 hours. For power computation, it 1s assumed that the variance of 


repairing time 1s 4.0. The appropriate hypothesis are 


Ho: p 2 40, 
Ha:p < 40. 


The experimenter decides to use a sample of n= 25 observations and @ =.05. He would 


like to know what the power of the t-test is for a range of possible population means. 


He decides to check means from 38 to 40 in .I increments. 


b. Inputs 
FIN ESOS: step 1: one sample one-sided t-test 
step 2: Wy = 40 
step 3: o = 2.0 
step 4: a@ =.05 


Step S2 sample sigesa 20 
step 6: maximum pt, =40., minimum ft, = 38, increment=.1 


c. Start program 


START PROGRAM: 
Do you want to test one sample t-test(y/n)? 


y 


Sample size ? 


7 


25 


Alpha-level ? 
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Figure 4.1a System flowchart for t-test power. 
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Figure 4.1b System flowchart for t-test power(cont’d.) 
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Figure 4.1c System flowchart for t-test power(cont’d.) 
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Figure 4.1d System flowchart for t-test power(cont’d.) 
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O5 
Mu-zero value ? 
2 
40 
Standard deviation value ? 
? 
? 
Do you want to use one-sided test(y/n)? 
y 
Do you want to test Ho: mu < mu0 vs Ha: mu > muO (y/n)? 
n 
Do you want to test Ho: mu 2 mu0 vs Ha: mu < muO (y/n)? 
y 
Maximum mul value (the maximum value on X-axis)? 
? 
40 
Minimum mul value (the minimum value on X-axis)? 
(condition : mul must be less than muO ) 
? 
38 
Increment mul value ? 
5 


d. Output 


The screen output is as follows: 


39 


LL, power 


NN= 38:000 SPO WER ae. 7. 

NN= 38.100 POWER= 0.99849 
NN= 38.200 POWER= 0.99677 
NN= 38.300 POWER= _ 0.99346 
NN= 38.400 POWER= _ 0.98742 
NN= 38.500 POWER= _ 0.97706 
NN= 38.600 POWER= _ 0.96028 
NN= 38.700 POWER= 0.93465 
NN= 38.800 POWER= _ 0.89772 
NN= 38,900" POWER See0rs47 op 
INN= 39:000> POWER >= 05/5526 
NN= 39.100" POW EK — 0 Uoe8 
INN 00 ae Oy tee cecal 
NN= 39:3005 “POWER = 052205 
NN= 39.400 POWER= _ 0.42564 
NN= 39.500 POWER= 0.33355 
NN= 39.600 POWER= _ 0.25053 
NN= 39.700 POWER= _ 0.17996 
NN= 39.800 POWER= _ 0.12338 
NN= 39.900 POWER= _ 0.08062 
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0.99933 fo X en Kea X enh pen Xan Kk enn pm nnn eee —————— $eereer--- +2 ------ +—+ 


power | s | 
| . | 

| . | 

| | 

| . | 

| : | 

| | 

| i | 

| | 

Ono oo9/ + : + 
| | 

| : | 

| | 

| : | 

| | 

| a | 

| | 

| : | 

| il 

Me OLOG=O] © taten--<--<- eee aphasia fas fone ace a ieee ses Re 

38. 00 38.95 I, 39-90 
X-SCALE: ”-”= 0.237E-01 UNITS 


Y-SCALE: *|"= 0.153E-01 UNITS 


2. Example of One sample two-sided t-test 
a. Scenario 
The supervisor in a rocket propellant plant is concerned with the burning 
rate of a rocket propellent. Standards specify that the mean burning rate must be 40 
inches per second. For power computation, it is assumed that o* = 4.0. The 
appropriate hypotheses are 
Ho: pu = 40 
Ha: pn # 40 


4] 


The engineer decides to use a sample of n= 25 observations and @ =.05. He would like 
to know the power of the t-test is for a range of possible population means. He 
decides to check means from 38 to 42 in .I| increments. 
b. Inputs 
INPUTS: step 1: one sample two-sided t-test 
step 2: Hy = 40 
StGO! 5140 ae 
step 4: @ =.05 
step 5: sample size = 25 
step 6: maximum ft; =42, minimum ff, = 38, increment=.1 
c. Start program 
START PROGRAM: 
Do you want to test one sample t-test(y/n)? 
y 
Sample size ? 
: 
25 
Alpha-level ? 
? 
OS 
Mu-zero value ? 
? 
40 
Standard deviation value ? 
? 
2 
Do you want to use one-sided test(y/n)? 
n 
Do you want to use two-sided test(y/n)? 
y 


a 


Maximum mul valve (the maximum value on X-axis)? 
9 


42 
Minimum mul value (the minimum value on X-axis)? 


(condition : mul must be less than mu-zero) 


42 


? 
38 
Increment mul value ? 
? 
nl 
d. Output 


The screen output is as follows: 


Ee 


NN= 
NN= 
NN = 
NN= 
NN= 
NN= 
NUN 
NN= 
NN= 
NN= 
NN= 
NN= 
NN= 
NN= 
NN= 


38.000 
38.100 
38.200 
38.300 
38.400 
38.500 
38.600 
38.700 
38.800 
38.900 
39.000 
39.100 
39200 
39.300 
39.400 
39.500 
39.600 
39.700 
39.800 
39.900 
40.000 
40.100 
40.200 
40.300 
40.400 
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POWER 


POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER = 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
CO Te 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
LOWERS 
POWER= 
POWER= 
POWER= 


U9 770 
99525 
OL SIONS 
0.98279 
0.96965 
0.94910 
0.91875 
0.87639 
082037 
Or ol09 
0.66944 
0.57882 
0.48379 
0.38975 
0.30192 
0.22459 
0.16066 
0.11146 
0.07708 
0.05692 
0.05028 
0.05687 
0.07698 
OFT! 
0.16045 


NN= 
NN = 
NN= 
NN= 
NN= 
NN= 


NN= 
NN= 


40.500 
40.600 
40.700 
40.800 
40.900 
41.000 
41.100 
41.200 
41.300 
41.400 
41.500 
41.600 
41.700 
41.800 
41.900 


4-5 


POWER= 
POW ER 
POWER= 
POWER = 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POWER= 
POV ie 
POWER= 
POWER= 
POWER= 


0.22433 
0.30162 
0.38942 
0.48344 
0.57848 
0.66912 
0.75081 
0.82034 
0.87621 
0.91862 
0.9490 I 
OF O99 
OS 27 9 
0.99070 
0.99524 


e297 70 Satie te +--------- +--------- +<—----—--- fase ER 


power | x x oe | 
| * . | 
| . * | 
| * i | 
| * f | 
| | 
| : : | 
| | 
| : | 
0. 52398 + 12 
| | 
| | 
| : . | 
| | 
| é * | 
| * ‘ | 
| | 
| . . | 
| kok tek | 
meow) Se-0) .+=+--ssas>-— Se So co Se ta SSS ———-=--== t+ 

38. 00 oa fH, 41.90 

X-SCALE: ”-”= 0.487E-01 UNITS 
Y-SCALE: ”|”= 0.158E-01 UNITS 
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V. DESCRIPTION OF THE POWER PROGRAM 


A. PROCEDURE OVERVIEW 

The interactive program included in this thesis is written in FORTRAN 77. It 
was written for use on an IBM370 from IBM3278 terminal. It is an interactive 
program. The command “NONCENT’ 1s all that is required to start this program. The 
executive file name imitializes the virtual machine environment and asks the user 
questions about the choice of program and compilation requirements. It then activates 
the selected program. Both ANOVA and TTEST are interactive programs which will 
be discussed in detail in later chapters. The output from the programs 1s presented on 


the terminal screen. 


£2. INSTRUCTIONS FOR PROGRAM ACCESS 

To start the program, make sure you have loaded (NONCEN, ANOVA, and 
TTEST) on your disk. In CMS (operating system mode), type “NONCEN’ and you will 
See. 

Please provide the FILENAME for your VS FORTRAN program. 
Now type the program name you want, for example “TTEST’, the response 1s: 
Do you need to compile your program ? (y/n) 

If you want to run, type “Y’ and your program will be loaded. 

A detailed description of the ANOVA (Analysis of Variance Test ) 1s given in 
Chapter III and of the t-test in Chapter IV. 

Following the screen output from the selected test, the user will be asked some 
questions about the output: 
Do you wish to BROWSE your output? (Y) 
n 
Print your output file? (Y) 
n 
Do you wish to XEDIT the program file? (Y/N) 
n 
Do you Wish to run the program again? (Y) 
n 


Then return to CMS mode. 
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VI. SUMMARY AND CONCLUSIONS 


Power considerations are useful in the design and assessment of statistical tests. 
The calculation of power usually involves noncentral distributions for which tables of 
probability are not available. A search was made for algorithms approximating the 
noncentral t, F and y 2 distributions. Algorithms giving sufficient accuracy and 
making efficient use of computer resources have been implemented in this thesis. 
Listings of the FORTRAN code for these implementations are included. 

An interactive program to compute and display power curves for several t-test 
and F-test situations has been developed. This program is user friendly and is 
described in this thesis. A listing of the program is provided. This program should be 


useful to researchers, experiment designers and statisticians. 
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AL Z END: 
AN INVERSE CENTRAL F APPROXIMATION 


Let P be the central F complementary CDF, and let Q be the 
complementary standard normal CDF [ Ref. 4:p. 947]. 


Define Vp by Q(Y_ =P, and Fy by PCF, | Vy Vo )=p. 


Then Fo ~ exp(2w) 


where w = wl-(w2)(w3) 

wl = (yp V (hth) )/h 

we = (1/(v 51 ee 

w3 = (h+(5/6)+2/(3h)) 

ho = 2/(C1/ (Vi eee ey 

h=(y, © -3)/6 

The program we used in the computation of critical values of F-tests is shown 

below. 


HHH HIKE HIM HH HEI HE HIE IEE HEHEHE HEHEHE EEE IESE EFE HEHE FEFEHE HE FEHE HEHE HEHEHE FE FEFEIEDE HE FE FEHE FEE PEDHEDEDEFHE HE FEE HEE FEE FE 


* % 
* APPROXIMATION TO THE F- INVERSE DISTRIBUTION. * 
* CALCULATE F-INVERSE GIVEN DOF1,DOF2,ALPHA. * 
x Q( F-ALPHA |DOF1,DOF2 )= ALPHA. * 
* * 
* * 


HHI III IK I IEE IE EHH IEE HIE IEICE EEE HIE SEE FEE FEES IEE FOIE IESE JE HEE IE IE EH GEIS FEE IIE EIEIE 
REAL ALPHA ,NUMER ,DENUM,LAM,xXP ,YP,FINVER,A>B 
DATA NF1,NF2,ALPHA/2 >4, .01/ 
DATA CO,C1,C2/2.515517,0.802853,0.010328/ 
DATA D1,D2,D3/1.432788,0.189269,0.001308/ 
B=NF1/2. 
A=NF2/2. 


CALCULATE NORMAL QUANTILE. 


T=SQRT(ALOG(1./( ALPHA®*2. ))) 
NUMER=CO+C1*T+C 2*Tx2 
DENUM=1.+D1*T+D2*T*2.+D3*T*x3. 
XP =T -( NUMER/DENUM ) 


OOO 


CALCULATE THE INVERSE F-DISTRIBUTION. 


LAM=( XP9€2-3. 0/6. 

W2= (1./(02.*B-1. ))-€1./(2.*A-1. )) 

WG= (1./02.%*A-1. ))401./02.*B-1. )) 

H=2 ° JWG 2636 ( -l ° } 

W1=(LAM#(5./6. )-(2./(3.*H))) 

W3=( XP%*(H+LAM )**.5)/H 

W=N3-WN2*W1 

FINVER=EXP( 2*W ) 

WRITE(3,100) NF1,NF2,ALPHA,FINVER 

100 eee ae NF1l=',1I3,‘'NF2=',I3,‘ALPHA=' »F5.4, ‘FINVERSE=',F6.3) 

0 

END 


aOoO 


The following table shows some values obtained using this approximation, 
together with exact values of critical values of the F distribution. It can be seen in 
Table § that this approximation provides adequate accuracy for the computation of the 
critical value of the F-tests. 
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TABLE 8 
mime wt LION ANID EXACT CRIMICAL VALUES OF F-TESTS 


approx 


4.150 
B25 
2.648 
eo 
1.668 
4.046 
O72 
2.804 
3) 
1.607 
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APPENDIX B 
THE INVERSE CENTRAL T APPROXIMATION 


We can use the asymptotic expansion for the Inverse t-CDF, as 
follows. Let Q be the complementary CDF of standard normal, and let P 
be the central t CDF [Ref. 4:p. 949]. 

ci P(t | v ) = 1-2p and (Xp ) =p 
then 

ty~ Xp * (91 (Xp )/ v)t (92 (xp) / ve )#(93 (Xp) / VP HC 94 OK) / VE” 
where” (97) (vm aco exe 

go (x)=(1/96)(5x°_+16x>_+3x), 
93 (x)=(1/384)( 3x7 #19x? +17x° = 5 ean) 
gq (x)=(1/92160)(79x? +776x’ +1482x° -1920x> -945x). 

The program we used in the computation of the critical value of the 
t-test is shown below. 

HHH HHI HEHEHE HEH HEH HEH HEHEHE IE HEE HE HE HEHEHE HE IEE HE HE HEHEHE EE EE EE HE EEE HE IEE HE HE EE EE HE HE IEE NE HE IE HE EIE IE HEHE HEN 
APPROXIMATION TO THE T- INVERSE DISTRIBUTION. x 
CALCULATE T-INVERSE GIVEN DOF ALPHA. * 
Q( T-ALPHA [| DOF )= ALPHA. x 

* 


HE HEHEHE HEE FETE HE IEE FE IE IE IE DE IEE ESE SE FE IESE FE FE FE FE IE FE IESE TE IE FEE HE HE EDE EEE IEE EME IEDEHEEDEDE EFE ESE HE FEE EE FE FEE HEHE FEE 


set initial conditions. 


QOD Kk KK KK XK 


REAL ALPHA ,NUMER ,DENUM,LAM,XP,G1,G62,G63,X,TINVER 
INTEGER N 

DATA N,ALPHA/9,.025/ 

DATA CO,C1,C2/2.515517,0.802853 ,0.010328/ 

DATA 01,02,03/1.432788,0.189269,0.001308/ 


CALCULATE NORMAL QUANTILE. 


Qa 


ALPHA=.05 
T=SQRTCALOG( 1. /CALPHA¥*2. )))} 
NUMER=CGO+C1*T +C2*xT x2 

DENUM=1. +D1*¥T4+D2*T**2.4+D3*%T*3. 
XP=T -(NUMER/DENUM ) 


CALCULATE THE INVERSE T-DISTRIBUTION. 


Nan 


X=XP , 

G1l=. 25*(X¥*3 +X ) 

G2=(1.796. J*¥(5. XXX*X54+16 . XXX*3 43. ¥X) 

G3=(1./384. )*(3 . XX*H7419. XX¥XSZ417 . XXHXSZ-15. ¥€X ) 

G4=(1./92160. )*0 79. XX*H94776 , HXHH741482 . XXXX5-1592 . RXKHS-9G5 , ¥X ) 

TP=X+(GI/NJ+0G2/0 NN ) 4+063/70N%*3 ) 4164/1 N¥*G ) ) 

TINVER=TP 

WRITE(3,100} N,ALPHA,TINVER 

100 o> DOF=' ,1354X, ‘ALPHA=' ,F5.3,6X,'TINVERSE='>F9.5) 

TO 

END 


The following table shows some values obtained using this 
approximation together with exact values of critical values of the t 
distribution. It can be seen in Table 9 that this approximation 
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—— 


provides adequate accuracy for computing the critical regions of the 
tc cests. 


WABEE 9 
APPROXIMATION AND EXACT CRITICAL VALUES OF T-TESTS 


DOF @ ( Type I error ) 


















































exact approx exact approx 

Z Z.919 Zee 6.964 6.475 
= Zao ZZ. 3.746 S27 01 
6 1.943 1.942 3.142 Se 3.0 
8 i oS3 1 ese Zee Zeros 
10 eo Z roZS 2.764 Ze Od 
eZ ie oS LSS 2,081 Ze 
14 iow ego 2.624 2.624 
6 1. 746 1. 746 Zee S ZS 
18 1.734 1.734 Aree PARES SI 
20 i= 2S in 25 ZaoZS Va es Fate) 
25 1.708 7 OS 2.485 2.485 
So 1. 690 ae © 2.438 2.438 
45 i 690 eo © ZZ 2.413 
65 i6G9 669 ZS oo Zs oS 
85 im ooS eos Bao yl ZS fk 
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APPENDIX C 
AN APPROXIMATION OF THE NORMAL CDF 


Let Z be the standard normal density function and let P be the 
corresponding CDF [Ref. 4:p. 939]. Then 


P(X) = f Z(t)dt ~ 1-Z(x)(azttayt@+ag,3)+ © (x) 
where t=1/(1+px) 
| © (x)| < .000001 


p =. 33267, and a, = . 4361836 
ay = — .1201676 
az = . 9372980 


D2 


APPENDIX D 
THE COMPUTATION OF NONCENTRAL F DISTRIBUTION 


HEHEHE HE HEHE IEFE HEHEHE HE IENE 36 HEHE HEHEHE «= INFORMATION 363636363636 3636 36363636 26 HE HEI HE HEHEHE HEHEHE HE HEHEHE HEHEHE HEHEHE 


THE OBJECTIVE OF THIS PROGRAM IS TO CALCULATE NON-CENTRAL F 
DISTRIBUTION (CDF ,F-INVERSE }. 
THE CALCULATION METHOD IS NORMAL APPROXIMATION. 


*K K KK OK XK 
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FVAL : X-VALUE DIMENSION 
FCDF : NON~CENTRAL F CDF DIMENSION. 
FINVER : GIVEN X VALUE, FIND INVERSE VALUE. 
FLAM : NON-CENTRAL PARAMETER 
FX : TO FIND INVERSE, GIVEN X VALUE. 
FSTAR : X VALUE VARIABLE. 
NF1l : DEGREE OF FREEDOM 1. 
NF2 : DEGREE OF FREEDOM 2. 
CDF : NORMAL DISTRIBUTION FUNCTION. 
P( 2 <= X) 
CONST : 1 7 SQRT ( 2 * PI ) IN NORMAL DISTRIBUTION FUNCTION. 
B : THE VALUE OF ‘X' IN NORMAL DISTRIBUTION. 


KK OK OK KK OK OK OK OK OK OK OK OK OK OK OK OK OK OK «OK 


FEIEIEH FE IEIEHE HE HEE HEPC HE HEIERE IEE HEHE HE HEFE IE HE HE IEE HE FE IE HEHEHE HEHE HEHE HEHE EIEDE HEHE HEHEHE HEHE HEIE IE HE IEIE IE IEE DE SE FE IE HEIEIE IE HEHE 
SET INITIAL CONDITIONS. 


REAL FVAL(100),FCDF({100),FINVER,FLAM,FMAX,FX,FINIAL 
DATA NF1,NF2,FLAM,FMAX/2>4,75.0,90.7 
DATA A1,A2,A3,P/.4361836,-.1201676,.9372980, .33267/ 
DATA FSTAR/25.128/7 

WRITE(3,2¢) NF1,NF2,FLAM 

WRITE(3,21) 

WRITE(3,22) 

WRITE(3,23) 

CONST = 1.0 / SQRT (2.0 * 3.1415927) 


CALCULATE THE NORMAL APPROXIMATION . 


C1l=(NF1*¥FSTAR )/(NF1+4FLAM ) 
C2=1-(2.7(09.*NF2 )) 
C3=(2.%*€NF1+2.*FLAM) 709. *(NF1+FLAM )**2. ) 
CC1l=(C1l¥*(1./3. ))*C2-(1-C3 ) 
CCO2=((C34+02.7(9. XNF2 ) (C1 02.73. ) 9.5) 
FX=CC1/CC2 
CDF = 0.0 
IF(FX .£€Q.0) THEN 
CDF = 0.5 
ELSE IF(FX.GT.O) THEN 
T=1.7(1.+PFX ) 
XSQU=( FX¥FX )*.5 
ZPDF=CONSTXEXP( -XSQU ) 
CDF=1-ZPDF*( ALT tAZ*THHZ+ASHT 23 ) 
ELSE 
FX=-FX 
T=1./(1. +PX¥FX ) 
XSQU=( FX¥FX J*.5 
ZPDF =CONST *EXP( -xXSQU ) 
CDF =ZPDF¥( AL¥T AZT H*2 4+ AZ KT HHS ) 
END IF 
FVAL(FSTAR)=FSTAR 
FCDF(FSTAR )=CDF 
WRITE(3,250 ) FSTAR,,FCDF(FSTAR) 
250 FORMAT('0O',F7.3,8X,F8.5) 
21 FORMAT('O',40('_')) 
Ze FORMAT( 'O' 52%, '"X"-VALUE ' ,8X%,'F(Z <= X)') 
2s FORMAT('0',40('_')) 
24 FORMAT('O','DOF1=',I2,3X,'DOF2=' ,I2,3X, 'NON-CENTRAL PARAMETER= 
&',F8.49%) 
STOP 
END 
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APPENDIX E 
THE COMPUTATION OF NONCENTRAL T DISTRIBUTION 


FEE IE EH HEHE HE IEEE IE EME IEEE HMIHENHM «INFORMATION 636363636363 IEE IIE IEE IEE IE ETE IE HE HEHE IESE IEEE EH HEHE 


THE OBJECTIVE OF THIS PROGRAM IS TO CALCULATE NON-CENTRAL T 
DISTRIBUTIONS ( CDFs) 

THE CALCULATION METHOD IS NORMAL APPROXIMATION. 

(REF. SEVERO & ZELEN'S (1960) ) 


KX KK KK XK 


HEHEHE IE HH ICH HH HHH IHHAHMHHHHHM VARIABLE DEFINITION ©3626 2636363696565 3 HH HHH III HHH 


TVAL : X-VALUE DIMENSION 
TCDF : NON-CENTRAL T COF DIMENSION. 
TLAM : NON-CENTRALITY PARAMETER 
TX : TO FIND INVERSE, GIVEN X VALUE. 
TSTAR : X VALUE VARIABLE. 
NF : DEGREE OF FREEDOM 
CDF : NORMAL DISTRIBUTION FUNCTION. 
P ( Z <= ®) 
CONST : 1 / SQRT ( 2 * PI J} IN NORMAL DISTRIBUTION FUNCTION. 
B > THE VALUE OF ‘X' IN NORMAL DISTRIBUTION. 


* KK KK KK KK KOK OK OK OK 


FEI HE HEHEHE HE HE HEHE EEE HEE EE HE HEE IE IE HEE EFEDE JE IEE HEHE FE HE HEHE EHE HE EEE HEHE HEHE HE HEHE HE FE HEHE FE HE FEE HEHE EE HEHE HE HEHE HE HEHE 


SET INITIAL CONDITIONS. 


000 KK KK KKK KK KK KK OK OK OK OK OK OK OK 


REAL TVAL(100),TCDOF(100),TINVER,TLAM,TMAX»>TX>TINIAL 
DATA NF ,TLAM/24,33.D078/ 
DATA A1,A2,A35P/.4361836 »-.12016765.9372980,.33267/ 
DATA TSTAR/24.4949/ 

WRITE(3,24¢) NF»TLAM 

WRITE(3,21) 

WRITE(3,22) 

WRITE( 3,23) 

CONST = 1.0 / SQRT (2.0 * 3.1415927) 


CALCULATE THE NORMAL APPROXIMATION . 
C1l=(1.-€1.704. NF) JXTSTAR-TLAM 


C2=SQRT(1. +0 TSTAR**2/02. *NF ))) 
TX=C1/C2 


000 


CDF = 0.0 
IF(TX .EQ.0)} THEN 
COF = 0.5 
ELSE IF(TX.GT.0) THEN 
T=1./01. +P*TX ) 
RSQU=(TX¥TX )*.5 
ZPDF =CONST EXP ( -XSQU } 
CDF=1-ZPDFX( A1L¥T tAZXT¥X2Z4AS%T *X3 ) 
ELSE 
TX=-TX 
T=1./01.+P*TX) 
RSQU=( TX¥TX )*.5 
ZPDF =CONST*EXP( -XSQU } 
CDOF=ZPDF( ALXT 4A2*T 2 4A3%T XS J 
END IF 
TVAL(TSTAR J=TSTAR 
TCDF( TSTAR )=COF 
WRITE(3,250) TSTAR,,TCDF(TSTAR) 
250 FORMAT('O',F7.3,8X,F8.5) 
21 FORMAT(‘0',40('_')) 
22 FORMAT( ‘0° ,2X,'"X"-VALUE ' »8X»,‘F(Z <= X)') 
23 FORMAT('O' ,4D('_')) 
24 FORMAT( 'D‘, ‘DOF=',I2,3X, 'NON-CENTRAL PARAMETER= 
&* ,F8.4) 
STOP 
ENO 
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APPENDIX F 
THE POWER OF ANOVA PROGRAM LIST 
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POWER ALGORITHM ( ANALYSIS OF VARIANCE }) 


DIRECTED BY : PROFESSOR DONALD.R.BARR 


WRITTEN BY : HUR, SEONG PIL JULY 1986 
DEPARTMENT OF O.R. NAVAL POSTGRADUATE SCHOOL 


* OK OK KK KOK OK KOK OK OK 


PROGRAM USES THE POWER PLOT TO TEST FOR ANALYSIS OF VARIANCE. 


(A) ONE-WAY ANOVA 
(1) SAMPLE SIZE VS POWER 
(2) NUMBER OF TREATMENT VS POWER 
(3) ALPHA-LEVEL VS POWER 
(4) NONCENTRALITY PARAMETER VS POWER 


(B) MULTI-WAY ANOVA 

(1) NUMBER OF REPLICATES VS POWER 

(3) ALPHA-LEVEL VS POWER 

(4) NONCENTRALITY PARAMETER VS POWER 
NOTE : subroutine PLOTT is NONIMSL subroutine library. 
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HHH HHI HHH HEHE HEHEHE IE HIE HEHEHE IEEE HE HEHEHE HEHE YEE EFL FE HEHE IEF FEE HEHEHE HEHEHE HEHEHE HEHE HEHEHE HEI 
REAL NN( 200) ,POWER( 200 ) ,NON,FINV ,NONCEN ,FPOW,FX,ALPHA,CDF ,CONST 
REAL MAXDEL »,MINDEL ,INCDEL ,STDS,PICON ,MINSD »MAXSD ,INCSD »MAXALP 
REAL MINALP ,INCRAL 
INTEGER MF( 20 ) »MAXN>,MINN>INCRN »>MAXK » MINK > INCRK 5 JJ 9K 9N 
CHARACTER*¥1 ANS 


START ONE-WAY OR MULTI-WAY ANALYSIS OF VARIANCE. 


Vekekekeke) 


00 PRINT*,'DO YOU WANT TO USE ONE-WAY ANOVA (Y/N)?' 
READ(5,7) ANS 
7 FORMAT(A1) 
IF((ANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 1000 
PRINT*¥,'DO YOU WANT TO PLOT N( # OF OBSERVATION PER TREATMENT ) 
& VS POWER ( Y/N)?' 
READ(5,7) ANS 
IF((ANS.EQ. ‘'Y').OR.(ANS.EQ. ‘Y')) GO TO 10 
ara eee YOU WANT TO PLOT K ( NUMBER OF TREATMENT ) VS POWER 
& ( Y/N)?" 
READ(5,7) ANS 
IF((CANS.EQ. ‘Y').OR.(CANS.EQ.‘Y')) GO TO 20 
PRINT*,'DO YOU WANT TO PLOT ALPHA-LEVEL VS POWER (Y/N)?' 
READ(5,7) ANS 
IF((ANS.EQ.‘'Y').OR.(ANS.EQ.‘Y')) GO TO 30 
PRINT*,'DO YOU WANT TO PLOT NONCENTRALITY VALUE VS POWER (Y/N)?' 
READ(5,7) ANS 
IF((ANS.EQ.'Y').OR.(ANS.EQ.‘'Y')) GO TO 40 
PRINT*¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ(557) ANS 
IFC CANS.EQ.'Y').OR.(ANS.EQ.'Y')) GO TO 500 
GO TO 999 


ANALYZE NUMBER OF SAMPLE SIZE PER TREATMENT VS POWER IN ONE-WAY 
ANALYSIS OF VARIANCE. 


FANNNANNO 


0 PRINT*,"NUMBER OF K( # OF TREATMENT )?' 
READ (5,*) K 
PRINT*, 'ALPHA- LEVEL? 
READ (5,*) ALPHA 
PRINT*,' STANDARDIZED SQUARED DISTANCE VALUE ?' 
READ (5,%*) STDS 


o> 


PRINT*,' MAXIMUM N VALUE ( THE MAXIMUM VALUE ON X-AXIS })?' 
READ (5,*) MAXN 
PRINT*, ‘MINIMUM N VALUE ( THE MINIMUM VALUE ON X-AXIS )? 
& (CONDITION :N MUST BE MORE THAN 2 )}' 
READ (5,*) MINN 
PRINT¥®,*INCREMENT N VALUE 7° 
READ (5,%) INCRN 
WRITE(6,95) 
95 FORMAT(1X,'DOF1 DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER‘) 
NF1=K-1 
JJ=1 
DO 110 N=MINN,MAXN, INCRN 
NF2=K¥(N-1) 


COMPUTE CENTRAL F-INVERSE GIVEN DOF1,DOF2 AND ALPHA-LEVEL. 


FINV=FINVER (NF1,NF2,ALPHA) 


COMPUTE NONCENTRALITY PARAMETER AS SAMPLE SIZE CHANGE. 


NONCEN=N¥STDS 


COMPUTE POWER USING NONCENTRAL F-DISTRIBUTION (CDF) GIVEN DOF1, 
DOF2,CENTRAL F-INVERSE AND NONCENTRALITY VALUE. 


MOMMONO ANNNN ANMINNN 


POW=FPO{ ALPHA ,NF1 »NF2»>F INV »NONCEN ) 
WRITE(6,15) NF1,NF2,ALPHA,FINV,NONCEN ,POW 
15 FORMAT( 'O',I5,3X%,I5,2X,2F10.551X,F12.251X,F10.5) 
NN JJ )=N 
POWER( JJ )=POW 
JJ=JJ+1 
110 CONTINUE 
JJ=JJ-1 
PRINT, ‘NN POWER ' 
DO 76 N=1>,JJ 
WRITE(6,75) NN(N),POWERIN ) 
75 FORMAT('O', ‘NN=',F8.3,4X, 'POWER=' ,F10.5) 
76 CONTINUE 


C 
C 
C PLOT NUMBER OF SAMPLE SIZE PER TREATMENT VS POWER IN ONE-WAY 
C ANALYSIS OF VARIANCE USING LIBRARY SUBROUTINE. 
C 
C 
CALL PLOTT{(NN,POWER,JJ,0) 
PRINT¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,7) ANS 
IF(CANS.EQ.‘Y').OR.CANS.EQ.'Y')) GO TO 500 
GOP 997 
C 
C 
C ANALYZE NUMBER OF TREATMENT VS POWER IN ONE-WAY ANALYSIS OF VARIANCE 
C 
C 
20 PRINT*,'# OF OBSERVATIONS PER TREATMENT (N= ) ?' 


READ (5,%*) N 
PRINT®, ‘ALPHA LEVEL (REAL VALUE) ?° 
READ (5,*%) ALPHA 
PRINT*,'WHAT STANDARDIZED SQUARED DISTANCE VALUE ?' 
READ (5,*) STDS 
PRINT, ‘MAXIMUM K (# OF TREATMENTS) RANGE ?'° 
READ (5,%) MAXK 
PRINT#, "MINIMUM K (# OF TREATMENTS) RANGE (K MUST BE MORE THAN 2)?' 
READ (5,%) MINK 
PRINT#,*°INCREMENT K VALUE ?' 
READ (5,%) INCRK 
WRITE(6,195) 

195 FORMAT(1X,'DOF1 DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER' ) 
JJ=1 
DO 111 K=MINK,MAXK,INCRK 
NF1=K-1 
NF2=K*(N-1 ) 
FINV=FINVER (NF1,NF2,ALPHA) 
NONCEN=N*STDS 
POW=FPO( ALPHA ,NF1,NF 2 »FINV ,NONCEN ) 
WRITE(6,115) NF1,NF2,ALPHA,FINV,NONCEN ,POW 

115 FORMAT('O',I5,3X,15,2X,2F10.5,1X,F12.251X»F10.5) 

NN( JJ )=K 
POWER{ JJ )=POW 
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JJ=JJ+1 

CONT INUE 

PRINT, ‘NN POWER ' 

DO 176 N=1>,JJ 

WRITE(6,175) NN(NJ,POWER(N ) 

FORMAT('O', ‘NN=',F8.354X, ‘POWER=',F10.5) 
CONT INUE 

CALL PLOTT(NN>,POWER >JJ>0) 

PRINT¥®,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,7) ANS 

IFC CANS.EQ.*Y').OR.CANS.EQ.'Y')) GO TO 500 
GO TO 999 


ANALYZE ALPHA-LEVEL VS POWER IN ONE-WAY ANALYSIS OF VARIANCE. 


PRINT*,'# OF OBSERVATIONS PER TREATMENT (N= ) ?' 
READ (5,*) N 

PRINT*,°'# OF TREATMENTS (K= )?' 

READ (5,%*) K 

PRINT*®,' STANDARDIZED SQUARED DISTANCE VALUE ?' 
READ (5,*) STOS 

PRINT*,'MAXIMUM ALPHA RANGE ?' 

READ (5,%) MAXALP 

PRINT*,'MINIMUM ALPHA RANGE 7?' 

READ (5,%*) MINALP 

PRINT*®, ‘INCREMENT ALPHA VALUE 7?' 

READ (5,*) INCRAL 
WRITE(6,295) 

FORMAT(1X,'DOF1L DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER’ ) 
JJ=1 
NF1=K-1 
NF2=Kx¥€N-1 ) 

DO 211 ALPHA=MINALP »MAXALP ,INCRAL 

FINV=FINVER (€NF1,NF2,ALPHA ) 
NONCEN=N*¥STOS 
POW=FPO( ALPHA »NF1,NF2,FINV,NONCEN ) 
WRITE(6,215) NF1,NF2,ALPHA,FINV,NONCEN,POW 
FORMAT( '0' ,I5,3X%,I5,2X%,2F10.5,1X,F12.2,1X,F10.5) 
NN( JJ )=ALPHA 
POWER( JJ )=POW 
JJ=JJ+1 
CONTINUE 
JJ=JJ—-1 

PRINT, ‘NN POWER * 

DO 276 I =1,JJ 

WRITE(6,275) NN€JJ),POWNER( JJ) 

FORMAT('O', ‘NN=',F8.3,4X, 'POWER=',F10.5) 
CONTINUE 
CALL PLOTT(NN>POWER,JJ,0) 
PRINT*,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,7) ANS 

IF(¢ANS.EQ. 'Y¥Y').OR.(ANS.EQ.‘'Y')) GO TO 500 
GO TO 999 


ANALYZE NONCENTRALITY PARAMETER VS POWER IN ONE-WAY ANALYSIS OF 
VARIANCE 


PRINT*,‘'# OF OBSERVATIONS PER TREATMENT (N= ) 7?' 

READ (5,*) N 

PRINT, '# OF TREATMENTS( K=)7' 

READ (5,*) K 

PRINT*, ‘ALPHA-LEVEL 7?' 

READ (5,%) ALPHA 

PRINT#, "MAXIMUM STANDARDIZED SQUARED DISTANCE VALUE RANGE ?° 
READ (5,%) MAXsD 

PRINT#, ‘MINIMUM STANDARDIZED SQUARED DISTANCE VALUE RANGE 7?' 
READ .€5,%)} MINSD 

PRINT#, ‘INCREMENT STANDARDIZED SQUARED DISTANCE VALUE RANGE 7' 
READ (5,%) INCSD 

WRITE( 6 »495 } 

FORMAT(1X,'DOF1 DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER' ) 
JJ=1 

NF1=K-1 

NF2=K*€N-1) 

DO 411 STDS=MINSD ,MAXSD,INCSD 

FINV=FINVER (€NF1,NF2,ALPHA) 

NONCEN=N*STDOS 

POW=FPOt ALPHA,NF1,NF2,FINV»,NONCEN ) 

WRITE(6,415) NF1,NF2,ALPHA,FINV,NONCEN,POW 

FORMATE '0' 515,3X,I5,2X,2F10.5,1X,F12.2,1X%,F10.5) 

NN( JJ )=NONCEN 
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POWER( JJ }=POW 

JJ=JJ+1 

CONTINUE 

JJ=JJ-1 

DO 476 I =1,JJ 

WRITE(6,4¢75) NN(I),POWER(T) 

FORMAT('O', '"NN=',F8.354X, 'POWER=',F10.5) 
CONT INUE 

CALL PLOTT(NN,POWER>JJ>0) 

PRINT¥,'OO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,7) ANS 
IF((CANS.EQ.‘Y').OR.(ANS.EQ.'Y'}} GO TO 500 
GO TO 999 


ANALYZE MULTI-WAY ANALYSIS OF VARIANCE. 


PRINT*, 'HOW MANY WAYS ANOVA CLASSIFICATION?’ 
READ(5,*) M 

DO 600 I=1,M 

PRINT*,‘FOR',I,'-FACTOR LEVEL »HOW MANY LEVELS ?' 
READ(5,*) MF(T) 

CONTINUE 


ANALYZE MODEL WITH FACTOR-LEVEL IN MULTI-WAY ANOVA. 


PRINT*,'IS THERE ANY INTERACTION ?(Y/N)° 
READ(5,7) ANS 
IF((ANS.EQ.'Y').OR.CANS.EQ. '¥Y')) GO TO 3100 
ISUM=0 
ITOT=1 
DO 800 I=1,M 

ISUM=ISUM+4(MF(I )-1) 

ITOT=ITOTXMF(T ) 
CONTINUE 
IEROF=ITOT-ISUM-1 


COMPUTE MINIMUM NUMBER OF REPLICATE PER TREATMENT. 


MINR=MINIR(ITOT ,ISUM ) 

PRINT*, ‘AT LEAST MINIMUM NUMBER OF REPLICATE PER TREATMENT =',MINR 
PRINT*,'WHICH FACTOR ANALYSIS (IF YOU WANT 3 RD FACTOR>PRESS 3)?' 
READ(5,*) I 

NF1=MF(I J-1 

PICON=FLOATCITOT )/FLOAT(MF(T)) 

GO TO 2000 


ANALYZE MODEL WITH FACTOR-LEVEL AND 2-WAY INTERACTION IN 
MULTI-WAY ANOVA. 


PRINT*,'ARE THERE ONLY 2-WAY INTERACTION AND FACTOR-LEVEL?(Y/N)' 
READ(5,7) ANS 
IF((ANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 2150 
ISUM1=0 
ISUM2=0 
ITOT=1 
DO 850 I=1,>M 
ISUM1=ISUM1+(MF(TI )-1) 
ITOT=ITOT*MF(T ) 
CONT INUE 
DO 860 I=1,M-1 
DO 870 J=I+1,M 
ISUM2=ISUM2+(MF(I J-1 )¥(MF(J)-1) 
CONTINUE 
CONTINUE 
ISUM=ISUM1+ISUM2 
IERDF=ITOT-ISUM-1 
MINR=MINIR(CITOT ,ISUM ) 
PRINT*¥,'AT LEAST MINIMUM NUMBER OF REPLICATE PER TREATMENT =',MINR 
PRINT*,'DO YOU WANT TO USE ONLY FACTOR-LEVEL (Y/N)?' 
READ(5,7) ANS 
IF((ANS.EQ.'N').OR.(CANS.EQ.'N')} GO TO 880 
PRINT*,'WHICH FACTOR ANALYSIS (IF YOU WANT 3 RD FACTOR,THEN 


&PRESS 3)?' 


READ(5,*) I 

NFI=MF(TI)-1 

PICON=FLOAT(ITOT J/FLOAT(MF(T )) 
GO TO 2000 


38 


880 PRINT*, WHICH INTERACTION ANALYSIS?( IF YOU WANT 2 AND 3 WAY INTER 
&ACTION,THEN PRESS 2 3 ).' 
READ(5,*) I>5J 
NFI=(MF(I J-1LJ*(MF(J)-1) 
PICON=FLOAT(ITOT )/FLOAT(MF(TI J¥MF(U)) 
GO TO 2000 


ANALYZE MODEL WITH FACTOR-LEVEL AND 3-WAY INTERACTION IN 
MULTI-WAY ANOVA. 


NAON90900 


150 PRINT*,'°ARE THERE ONLY 3-WAY INTERACTION AND FACTOR-LEVEL ?7(Y/N}' 

READ(5,7} ANS 

IF((ANS.EQ.'N'}.OR.(CANS.EQ.'N')) GO TO 2200 

ISUM1=0 

ISUM3=0 

ITOT=1 

DO 871 I=1,M 

ISUM1=ISUM1+(MF(I)-1) 
ITOT=ITOT*MF(T) 

871 CONT INUE 

DO 872 I=1,M-2 

DO 873 J=I+1,M-1 
DO 874 K=J+1,M 
ISUM3=ISUM3+(MF(I J-1 }*(MF( J )-1 J*0MF(UK J-1) 

874 CONTINUE 
873 CONTINUE 
872 CONTINUE 

ISUM=ISUM1+ISUM3 

IERDF=ITOT-ISUM-1 

MINR=MINIR(ITOT ,ISUM ) 

PRINT*,'AT LEAST MINIMUM NUMBER OF REPLICATE PER TREATMENT =',MINR 

PRINT*¥,'DO YOU WANT TO ANALYSIS ONLY FACTOR-LEVEL ?(Y/N)' 

READ(5,7) ANS 

IF(( ANS. EQ.'N’).OR.(ANS.EQ.‘'N')) GO TO 875 

PRINT®,*WHICH FACTOR ANALYSIS (IF YOU WANT 3 RD FACTOR>PRESS 3)?' 

READ(5,*) I 

NF1=MF(I)-1 

PICON=FLOATC(ITOT )/FLOAT(MF(I)) 

GO TO 2000 
875 PRINT*,'DO YOU WANT TO ANALYSIS 3-WAY INTERACTION ?(Y/N)' 
READ(5,7) ANS 
IF((ANS.EQ.'N').OR.(ANS.EQ.'N'))} GO TO 500 
PRINT®,'WHICH INTERACTION ANALYSIS?( IF YOU WANT 2 AND 3 AND 4 INT 
&ERACTION ,THEN PRESS 2 3 4 ).' 
READ(5,*) I,J,K 
NF1=(MF(I }-1 )*(MF( J )-1 *(MF(K )-1) 
PICON=FLOAT(CITOT )/FLOAT(MF(I )*MF( J *MF(K )) 
GO TO 2000 


ANALYZE MODEL WITH FACTOR-LEVEL ,2-WAY INTERACTION IN AND 
3-WAY INTERACTION IN MULTI-WAY ANOVA. 


NOQ0N00 


200 PRINT¥,'ARE THERE 2-WAY INTERACTIONS AND 3-WAY INTERACTIONS AND 
&FACTER-LEVEL (Y/N)?! 
READ(5,7) ANS 
IF( (ANS.EQ.'N'}.OR.( ANS.EQ.'N')) GO TO 999 
ISUM1=0 
ISUM2=0 
ISUM3=0 
ITOT=1 
DO 950 I=1)M 
ISUML=ISUM1+(MF(I)-1) 
ITOT=ITOTXMF(TI ) 
950 CONTINUE 
DO 960 I=1,M-1 
DO 970 J=I+15M 
ISUM2=ISUM2+(MFCI )-1)*(MF(J)-1) 
970 CONTINUE 
960 CONTINUE 
DO 965 I=1,M-2 
DO 975 J=I+1,M-1 
DO 985 K=J+1,M 
ISUM3=ISUM3+(MF(I )-1 )%*(MF(U )-1 )%*(MF(K )-1) 
985 CONTINUE 
975 CONTINUE 
965 CONTINUE 
ISUM=ISUM1+ISUM2+ISUM3 
IERDF=ITOT-ISUM-1 
MINR=MINIR( ITOT ,ISUM) 
PRINT*,'AT LEAST MINIMUM NUMBER OF REPLICATE PER TREATMENT ='»,MINR 
PRINT¥,'DO YOU WANT TO ANALYSIS ONLY FACTOR-LEVEL (Y/N)?! 
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READ( 5,7) ANS 

IF( (ANS.EQ.'N')J.OR.¢(ANS.EQ.'N')) GO TO 990 

PRINT*, ‘WHICH FACTOR ANALYSIS (IF YOU WANT 3 RD FACTOR,PRESS 3)?' 
READ(5,*) I 

NF1=MF(I)-1 

PICON=FLOAT(ITOT )/FLOAT(MF(TI )) 

GO TO 2000 

PRINT*,'DO YOU WANT TO ANALYZE 2-WAY INTERACTION (Y/N)?! 
READ( 5,7) ANS 

IF((ANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 995 

PRINT*, ‘WHICH INTERACTION ANALYSIS 7( IF YOU WANT 2 AND 3 WAY INTE 
&8ACTION,THEN PRESS 2 3 )}. 

READ(5,*) I,J 

NF1l=(MF(I )-1LJ¥¢(MFOCJUJ-1) 

PICON=FLOAT(ITOT )/FLOAT(MF( I )¥MF( J) ) 

GO TO 2000 

PRINT*¥,'DO YOU WANT TO ANALYZE 3-WAY INTERACTION (Y/N)?! 
READ( 5,7) ANS 

IF( (ANS.EQ. 'N‘').OR.(ANS.EQ.'N'))} GO TO 500 

PRINT*, ‘WHICH INTERACTION ANALYSIS?( IF YOU WANT 2 AND 3 AND G INT 
&ERACTION ,THEN PRESS 2 3 @ ). 

READ(5,*) I» Jy K 

NFLl=(MFC(I )-1 (MFC J J=-1 eC MFOK )-1) 

PICON=FLOAT(ITOT )/FLOAT(MF( I )¥MF( J )*MFUK ) ) 

GO TO 2000 

PRINT*¥,'DO YOU WANT TO PLOT R( # OF REPLICATE PER CELL) 
&VS POWER ( Y/N)?' 

READ(5,7) ANS 

IF( (ANS.EQ.'Y').OR.¢(ANS.EQ.'Y')) GO TO 4100 

PRINT¥,'DO YOU WANT TO PLOT NONCENTRALITY VS POWER (Y/N)?' 
READ( 5,7) ANS 

IF( (ANS.EQ.'Y').OR.(ANS.EQ.'Y')) GO TO 4200 

PRINT*¥,'DO YOU WANT TO PLOT ALPHA-LEVEL VS POWER (Y/N)?' 
READ( 5,7) ANS 

IF(( ANS.EQ.'Y'}J.OR.(ANS.EQ.'Y')) GO TO 2300 

PRINT*¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 

READ (5,7) ANS 

IF( (CANS.EQ.'Y').OR.¢(ANS.EQ.'Y’)) GO TO 500 

GO TO 999 


ANALYZE NUMBER OF SAMPLE SIZE PER TREATMENT VS POWER IN MULTI-WAY 
ANALYSIS OF VARIANCE. 


PRINT*, ‘ALPHA LEVEL?' 

READ (5,*)} ALPHA 

PRINT*¥, ‘STANDARDIZED SQUARED DISTANCE VALUE ?’ 
READ (5,*) STDS 

PRINT*,'MAXIMUM R VALUE ( THE MAXIMUM VALUE ON X-AXIS )?' 
READ (5,*) MAXR 

PRINT*¥,° MINIMUM R VALUE ( MUST BE MORE THAN',MINR,‘' )?° 
READ (5,%*) MINR 

PRINT*¥,° INCREMENT R VALUE ?' 

READ (5,*) INCRR 

WRITE(6 ,4095) 

FORMAT(1X,'DOF1 DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER‘ ) 
JJ=1 

DO 4110 N=MINR,MAXR,INCRR 

NF 2=N*¥ITOT-ISUM-1 

FINV=FINVER (NF1,NF2,ALPHA ) 

NONCEN=N¥PICON*¥STDS 

POW=FPO( ALPHA ,NF1,NF2,FINV,NONCEN ) 

WRITE(6,4¢015) NF1,NF2,ALPHA,FINV,NONCEN,POWN 
FORMAT( '0O',I5,3X%,I5,2X%,2F10.5,1X,F12.2,1X%,F10.5) 
NN( JJ )=N 

POWER( JJ )=POW 

JJeJJe1 

CONTINUE 

JJ=JJ=-1 

PRINT, ‘NN POWER ' 

DO 4076 N=1,JJ 

WRITE(6;4¢075) NN(N),POWER(N) 

FORMAT('O', '‘NN=',F8.35,4X%, 'POWER=',F10.5) 
CONTINUE 

CALL PLOTT(NN,POWER,JJ,0) 

PRINT*,‘DO YOU WANT TO PLOT AGAIN (Y/N)?' 

READ (5,7) ANS 

IF((ANS.EQ.'Y').OR.(ANS.EQ.‘'Y')) GO TO 500 

GO TO 999 


ANALYZE NONCENTRALITY PARAMETER VS POWER IN MULTI-WAY 
ANALYSIS OF VARIANCE. 
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C 
4200 ee MANY NUMNER OF REPLICATE (MUST BE MORE THAN’ ,MINR,> 
&' )7 
READ(5,*) IRR 
NF2=IRR¥ITOT-ISUM-1 
PRINT*¥,'ALPHA-LEVEL ?' 
READ (55%) ALPHA 
PRINT*», "MAXIMUM STANDARDIZED SQUARED DISTANCE VALUE RANGE ?'° 
READ (5,*) MAXSD 
PRINT¥,'MINIMUM STANDARDIZED SQUARED DISTANCE VALUE RANGE ?' 
READ (5,*) MINSD 
PRINTx, ‘INCREMENT STANDARDIZED SQUARED DISTANCE VALUE RANGE ?' 
READ (5,*) INCSD 
WRITE! 65,4495 ) 
4495 FORMAT(1X,'DOF1 DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER' ) 
JJ=1 
DO 4911 STDS=MINSD »MAXSD,INCSD 
FINV=FINVER (NF1,NF2,ALPHA ) 
NONCEN=IRR*¥PICON*STDS 
POW=FPO( ALPHA ,NF1,NF2,FINV,NONCEN ) 
WRITE(6,4915) NF1,NF2,ALPHA,FINV,NONCEN,POW 
4915 FORMAT('O',I5,3X,15,2X%,2F10.5,1X,F12.2,1X,F10.5) 
NN( JJ J=NONCEN 
POWER( JJ )=POW 
JJ=JJ+1 
4911 CONTINUE 
JJ=JJ—-1 
DO 1476 I =1,JJ 
WRITE(6,14¢75) NN(1I),POWER(T) 
1475 FORMAT('O', 'NN=',F8.354X, 'POWER=',F10.5) 
1476 CONTINUE 
CALL PLOTT(NN,POWER,JJ,0 ) 
PRINT¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,7) ANS 
IF(( ANS.EQ.'Y').OR.(ANS.EQ.‘'Y¥Y')) GO TO 500 
GO TO 999 


C 

c 

c ANALYZE ALPHA-LEVEL PARAMETER VS POWER IN MULTI-WAY 
C ANALYSIS OF VARIANCE. 
c 
C 
2 


300 PRINT®,'HOW MANY NUMNER OF REPLICATE (MUST BE MORE THAN’ ,MINR> 
&')?' 


READ(5,*) IRR 
NF 2=IRR*¥I TOT-ISUM-1 
READ (5,*) K 
PRINT*, ‘STANDARDIZED SQUARED DISTANCE VALUE ?' 
READ (5,*) STDS 
PRINT¥, ‘MAXIMUM ALPHA RANGE ?'° 
READ (5,%) MAXALP 
PRINT*¥,'MINIMUM ALPHA RANGE 7?'° 
READ (5,%*) MINALP 
PRINT*¥, ‘INCREMENT ALPHA VALUE ?' 
READ (5,%*) INCRAL 
WRITE(6,1295) 
1295 FORMAT(1X,'DOF1 DOF2 ALPHA-LEVEL F-INVERSE NONCENTRAL POWER’ ) 
JJ=1 
NF1=K-1 
NF2=K¥(N-1 ) 
DO 1211 ALPHA=MINALP ,MAXALP ,INCRAL 
FINV=FINVER (NF1,NF2>,ALPHA ) 
NONCEN=IRR¥PICON¥STDS 
POW=FPO( ALPHA »NF1l,NF2,FINV ,NONCEN ) 
WRITE(6,1215) NF1l»NF2,ALPHA,FINV,NONCEN,POW 
1215 FORMAT('0' ,I5,3X,I5,2X%,2F10.5,1X,F1l2.25,1X,F10.5) 
NN( JJ )=ALPHA 
POWER( JJ )=POW 
JJ=JJ4#1 
1211 CONTINUE 
JJ=JSS=1 
PRINT#> "NN(JJJ)° ,MINK ,»MAXK , INCRK 
DO 1276 I =1,JJ 
WRITE(6,1275) NN(JJ),POWER( JJ) 
1275 FORMAT('O', 'NN=',F8.354X, 'POWER=',F10.5) 
1276 CONTINUE 
CALL PLOTT(NN,POWER,JJ,0) 
PRINT*,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,7) ANS 
IF((ANS.EQ.'Y').OR.(ANS.EQ. 'Y¥')) GO TO 500 


GO TO 999 
999 STOP 

END 
JOH HEI IE IE IE IE IEEE EE IEICE TEE IE IEICE IEIE IE HEE IETHE IEEE HE HEE IESE EDIE DEE HE FETE IE IESE FEE TEE EEE HE SESE FE ESE HE FEFDEDE FE DE HE FEE FE FE 
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APPROXIMATION TO THE F- INVERSE DISTRIBUTION. 


CALCULATE F-INVERSE GIVEN DOF1,DOF2,ALPHA. 
Q{ F-ALPHA |DOF1,DOF2 )= ALPHA. 


HK HHH HK HIGHT KH HE HEHEHE HE HEHEHE HEHEHE HEHE HEHEHE HEHEHE HE HEHEHE HEHEHE FE MEE HEHE FIEDE HE HE FEFE FEE EEE FE FEEFE HE FE HE_ HEHEHE 


FUNCTION FINVER(N1,N2Z,ALP ) 

REAL ALP ,NUMER ,DENUM,LAM,XxXP ,YP,A,8,FINI 
DATA CO,C1,C2/2.515517,0.802853 ,0.010328/ 
DATA D1,D2,D3/1.432788,0.189269,0.001308/ 
B=N1/2. 

A=N2/2. 


CALCULATE NORMAL QUANTILE. 


T=SQRT(ALOG(1./( ALP¥*2. ))) 
NUMER=CO+C1¥T +C2xT HZ 

DENUM=1 . +D1*¥T+D2*T XX 24+D3%T HHS 
XP =T -( NUMER/DENUM ) 


CALCULATE THE INVERSE F-DISTRIBUTION. 


LAM=( XP¥*2-3. 7.76. 
W2= (1.702. *B-1. J) )-( 
WNG= (1.712. #A-1. ) +f 
H=2.*WG**(-1. ) 
W1=( LAM#+(5./6. )-(2./(03.*H) )) 
W3=(XP*(H+LAM )¥*.5 )/H 
W=(W3-NZ*W1 )*2Z 

IF(W.LT.174¢) THEN 

FINI=EXP(W) 


E 
FIN1=1000 
ENDIF 
FINVER=FIN1 
RETURN 
END 


(2.*A-1.)) 


1.7 
1.72. *6-1, 3) 


rs 
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GIVEN DOF1,DOF2,FINV,NON-CENTRAL PARAMETER; 
THEN COMPUTE POWER. 


( I.E. CALCULATE NON-CENTRAL F-DISTRIBUTION ( C.D.F.} 


HHH HH HH HH HH HHH HHH HHH HHH HHH HE HH HH EH HH HH HHI IEE HEIDE HE HEIEFE HEHE HEHE HEIEHE HEHE HEHEHE 


FUNCTION FPO( ALP »N1,N2,FIN,NON) 
REAL NON 
INTEGER N1;,N2 


DATA A1,A2,A3,P/.4361836,-.1201676,. 9372980, .33267/ 


CONST = 1.0 7 SQRT (2.0 * 3.1415927) 


CALCULATE THE NORMAL APPROXIMATION . 


C1=(N1*FINJ/(N1+NON ) 

C2=1-(2./19.¥N2)) 
C3=(2.%*(N1+2. NON ) )/(9.% (N1+NON ) 2 ) 
CC1=(C1¥¥(1./3. ))*¥C2-(1-C3 ) 
CCO2=((C34(2.709. NZ) CCL HI 2.7/3. ) d**.5) 
FX=CC1/CC2 


CALCULATE THE NONCENTRAL F-DISTRIBUTION USING NORMAL CDF. 


CDF = 0.0 
IF(FX .EQ.0) THEN 
CDF = 0.5 
ELSE IF(FX.GT.0) THEN 
T=1./(1. +PXFX ) 
XSQU=( FX*FX )%*.5 
IF(XSQU.LT.174) THEN 
ZPDF =CONST*EXP ( -XSQU ) 


E 
ZPDF=0 

ENDIF 
ZPDF=CONST *EXP ( -XSQU ) 
CDF =1-ZPDFX¥( ALXT tA2HT HH 2+ ASHT 43 ) 
ELSE 
FX=-FX 
T=1./(1. tP¥FX) 
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XSQU=( FX*FX J*.5 
IF(XSQU.LT.174¢) THEN 
ZPDF=CONST*EXP( -XSQU ) 


CDF=ZPDFX( ALXT +AZ*T X¥2+AS%T HHS ) 
END IF 
FPO=1-CDF 
RETURN 
END 
HII IISA HIE IEE IESE FEED IEE HEHE FEIEE HEE EEE DES SE TEI TETE ESET TE TEE TFC TEE FOIE ILI 
« x 
¥ GIVEN ITOT AND ISUM , THEN CALCULATE MINIMUM REPLICATE. * 
” x 
HEHEHE HEHEHE HEHEHE IEE IEIEFEIE IE PEE HEIDE FEI IEDESE FEI IEDEDEFEIEDEESEIESE BE FE IEDE DEDEDE DED FETE DESE FEE DEE IEEE IIE IEE 
FUNCTION MINIR(ITOT ,ISUM) 
INTEGER IR», ITOT,ISUM,ICHECK ,MINIR 
IR=2 
10 ICHECK=IR*ITOT-ISUM-1 
IF(ICHECK.LE.O) THEN 
IR=IR+1 
GO TO 10 
ELSE 
MINIR=IR 
ENDIF 
RETURN 
END 
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APPENDIX G 


THE POWER OF T-TEST PROGRAM LIST 
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PROGRAM IS USED THE PRESENTATION 


POWER ALGORITHM ( T TEST ) 
DIRECTED BY : PROFESSOR DONALD.R. BARR 


WRITTEN BY : HUR, SEONG PIL JULY 1986 
DEPARTMENT OF O.R. NAVAL POSTGRADUATE SCHOOL 


(1) ONE-SAMPLE CASE TEST. 
- ONE-SIDED TEST ( HA : MU>MUO ) 
- ONE-SIDED TEST € HA : MU<MUO ) 
- TWO-SIDED TEST ( HA : MU MUO ) 
(2) TWO-SAMPLE CASE TEST. 
- ONE-SIDED TEST € HA : MU>MUO ) 
- ONE-SIDED TEST (¢ HA MU<MUO ) 
- TWO-SIDED TEST ( HA MU MUO 


) 
NOTE : Subroutine PLOTT is NONIMSL subroutine library. 


HEH HH HH HHH HIE HIEKE HEHE HEHEHE HEHEHE IE HEHEHE HEHE DE HEHE HEE HEHE HE FEI TE IE HEE HE FE IE FE DE FE IE HEHE DHE IE IE HEE IE SEE DEE HE IE SHEE IEE 


REAL NN( 200 ),POWER( 200 ),DF ,NONCDF ,TCDF1,TCDF2,TINVI 

REAL TINV,NONCEN,TX,ALPHA,CDF ,CONST ,»MAXMU1 ,MINMU1 ,INCMU1 
REAL MAXDEL ,MINDEL ,INCDEL,POW,POW1,POW2 ,MUZERO,STDV ,MUL] 
INTEGER JJ>K>N 

CHARACTER¥1 ANS 


INPUT USER OPTIONS. 


PRINT®,'DO YOU WANT TO TEST ONE SAMPLE T-TEST(Y/N)?' 
READ(5,10) ANS 

FORMAT(AL) 
IF((ANS.EQ.'N*).OR.(ANS.EQ.'N‘)) GO TO 160 
PRINT*¥,' SAMPLE SIZE ?' 

READ(5,*) N 

PRINT*,' ALPHA-LEVEL ?' 

READ(5,*) ALPHA 

PRINT*¥," MU-ZERO VALUE ?' 

READ(5,*) MUZERO 

PRINT*,' STANDARD DEVIATION VALUE ?' 
READ(5,*) STDV 


ONE-SAMPLE ONE-SIDED TEST, HA : MU > MUO. 


PRINT¥,'DO YOU WANT TO USE ONE-SIDED TEST (Y/N)?' 

READ(5,10) ANS 

IF((CANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 30 

PRINT*,*DO YOU WANT TO TEST HO:MU=MUO VS HA:MU>MUO (Y/N)?' 
READ(5,10) ANS 

IF(CANS.EQ.'N').OR.(CANS.EQ.'N')) GO TO 40 

PRINT," MAXIMUM MU1 VALUE ( THE MAXIMUM VALUE OF X-AXIS )7?' 
READ (55%) MAXMU1 

PRINT#%,° MINIMUM MU1 VALUE ( THE MINIMUM VALUE OF X-AXIS )? 


& (CONDITION :MU1 MUST BE MORE THAN MU-ZERO )' 


READ (55%) MINMU1 

PRINT#,* LARGE INCREMENT N VALUE ?' 
— INCMU1] 

N=N- 


GENERATE T-INVERSE VALUE ( CRITICAL REGION VALUE ). 


TINV=TINVER(N,ALPHA ) 
JJ=1 
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*K KK KK KK OK OK OK OK OK OK OK KOK OK OK OK OK OK 


NNMaAaAN 


‘elelekels 


50 


ANAAaAN 


2 


‘elekeleke) 


70 
60 


40 
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120 
110 


30 


DO 50 MU1L=MINMU1 »MAXMU1 , INCMU1] 


COMPUTE NON-CENTRALITY VALUE. 


DF=N 
NONCEN=( (MU1-MUZERO )*SQRT (DF +1) )/STDV 


GENERATE CDF OF NONCENTRAL T DISTRIBUTION. 


TCDF=NONCDF (N,>,NONCEN,TINV ) 

POW=1-TCDF 

NN( JJ )=MU1 

POWER( JJ )=PON 

JJ=JJ+1 

CONT INUE 

JJ=JJ-1 

PRINT,» ‘NN POWER ' 

DO 60 I=1,JJ 

WRITE(6,70) NN(I),POWER(T ) 

FORMAT('O', '‘NN=',F8.354X, '"POWER=',F10.5) 
CONT INUE 

CALL PLOTT(NN,POWER,JJ,0) 

PRINT¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,10) ANS 
IF((ANS.EQ.'Y').OR.(ANS.EQ.'Y¥')) GO TO 100 
GO TO 80 


ONE-SAMPLE ONE-SIDED TEST, HA : MU < MUO. 


PRINT¥,'DO YOU WANT TO TEST HO:MU=MU0O VS HA:MU<MUO (Y/N)7?' 
READ(5,10) ANS 

IF((ANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 80 

PRINT%,' MAXIMUM MU1 VALUE ( THE MAXIMUM VALUE OF X-AXIS })?' 
READ (5,%*) MAXMU1 

PRINT*,' MINIMUM MU1 VALUE ( THE MINIMUM VALUE OF X-AXIS }? 


& (CONDITION :MU1 MUST BE LESS THAN MU-ZERO ) ' 


READ (5,*) MINMU1 
PRINT¥,° INCREMENT MU1 VALUE ?' 
READ (5, JINCMU1 
N=N-1 
TINV=TINVER(N,ALPHA } 
TINV=-TINV 
JJ=1 
DO 90 MU1L=MINMU1 ,MAXMU1 , INCMU1 
DF=N 
NONCEN=( (MU1-MUZERO )¥SQRT (DF +4+1))/STDV 
TCDF=NONCDF(N,NONCEN, TINV ) 
POW=TCDF 
NN( JJ )=MU1 
POWNER( JJ )=POW 
JJ=JJ+1 
CONT INUE 
JJ=JJ-1 
PRINT, 'NN POWER‘ 
DO 110 I=1,JJ 
WRITE(6,120) NN(1I),PORER(T ) 
FORMAT('O', '‘NN=',F8.3,4¢X%, 'POWER=',F10.5) 
CONTINUE 
CALL PLOTT(NN,POWER,JJ,0) 
PRINT¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (55,10) ANS 
IF((ANS.EQ.'Y'}.OR.(ANS.EQ.'Y'})} GO TO 100 
GO TO 80 


ONE-SAMPLE TWO-SIDED TEST, HA : MU =/ MUO. 


PRINT¥,'DO YOU WANT TO USE TWO-SIDED TEST(Y/N)?' 

READ.5,10) ANS 

IF((ANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 80 

PRINT, ‘MAXIMUM MUL] VALUE ( THE MAXIMUM VALUE ON X-AXIS }7?' 
READ (5,%*) MAXMU1 

PRINT¥,* MINIMUM MU1 VALUE ( THE MINIMUM VALUE ON X-AXIS )? 


& (CONDITION :MU1 MUST BE LESS THAN MU-ZERO }' 


READ (5,%) MINMU1 

PRINT¥,° INCREMENT MU1 VALUE ?' 
READ (5,% )INCMU1 

N=N-1 
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130 


150 
140 
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C 


190 


210 
200 


RFONNAN 


ALPHA=ALPHA/Z2. 
cde en ace 
=1 
ae 130 MUL=MINMU1 ,>MAXMU1 , INCMU1 
DF=N 
NONCEN=( (MU1-MUZERO )*SQRT(DF+1. ))/STDV 
TCDF 1=NONCDF (N,NONCEN, TINV ) 
POW1=1.-TCDF1 
TINV1=-TINV 
TCDF 2=NONCDF(N,NONCEN >TINV1 ) 
POW2=TCDF2 
POW=POW1+POWZ 
NN( JJ) =MU1 
POWER( JJ )=POW 
JJ=JJ+1 
CONTINUE 
JJ=JJ~-1 
PRINT, "NN POWER ' 
DO 140 I=1,JJ 
WRITE(6,150) NN(I),POWER(T) 
FORMAT('O', 'NN=',F8.3,4X, 'POWER=',F10.5) 
CONTINUE 
CALL PLOTT(NN,POWER,JJ,0) 
PRINT¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,10) ANS 
IF((CANS.EQ.'Y').OR.(ANS.EQ.‘'Y')) GO TO 100 
GO TO 80 


TWO-SAMPLE ONE-SIDED TEST, HA : MUX > MUY. 


PRINT*¥,'DO YOU WANT TO TEST TWO SAMPLE T-TEST(Y/N)?' 
READ(5,10) ANS 

IF((CANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 80 

PRINT®,* SAMPLE SIZE FROM POPULATION 1 ?' 

READ(5,*) N 

PRINT¥,' SAMPLE SIZE FROM POPULATION 2 ?' 

READ(5,%*) M 

PRINT*,' ALPHA-LEVEL 7?' 

READ(5,*) ALPHA 

PRINT¥,'DO YOU WANT TO USE ONE-SIDED TEST (Y/N)?' 
READ(5,10) ANS 

IF((ANS.EQ.'N' ).OR.(ANS.EQ.'N')) GO TO 170 
PRINT¥,'DO YOU WANT TO TEST HO: MUX=MUY VS HA: MUX>MUY (Y/N)?' 
READ(5,10) ANS 

IF(CANS.EQ.'N').OR.(ANS.EQ.‘'N')) GO TO 180 

PRINT¥,' MAXIMUM NONCENTRALITY VALUE ? ' 

READ (5,%) MAXDEL 

PRINT*¥,' MINIMUM NONCENTRALITY VALUE? ' 


& (CONDITION :MU1 MUST BE MORE THAN MU-ZERO )7' 


READ (5,%) MINDEL™ 

PRINT*®,' INCREMENT NONCENTRALITY VALUE?' 
READ (5,*jJINCDEL 

N1=M+N~-2 

TINV=TINVER(N1,ALPHA ) 

JJ=1 

DO 190 NONCEN=MINDEL »MAXDEL »INCDEL 
TCDF=NONCDF(N1,NONCEN,TINV ) 

POW=1-TCDF 

NN( JJ )=NONCEN 

POWER( JJ )=POW 

JJsJJ+1 

CONTINUE 

JJ=JJ-1 

PRINT, 'NN POWER * 

DO 200 I=1,JJ 

WRITE(6,210) NN(I),POWER(T) 

FORMAT('O', '‘NN=',F8.3,4X,'POWER=',F10.5) 
CONTINUE 

CALL PLOTT(NN,POWER,JJ,0) 

PRINT®,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,10) ANS 

IF((CANS.EQ.'Y' ).OR.(ANS.EQ.'Y')) GO TO 100 
GO TO 80 


TWO-SAMPLE ONE-SIDED TEST, HA : MUX < MUY. 


PRINT¥,'DO YOU WANT TO TEST HO:MUX=MUY VS HA:MUX<MUY (Y/N)?° 
READ(5,10) ANS 

IF((ANS.EQ.'N').OR.(ANS.EQ.'N')) GO TO 80 

PRINT*,' MAXIMUM NONCENTRALITY VALUE 7? ' 

READ (5,%) MAXDEL 
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PRINT*,* MINIMUM NONCENTRALITY VALUE? ' 
REAO (5,x*) MINDEL 

PRINT¥®,* INCREMENT NONCENTRALITY VALUE?‘ 
REAO (5, )JINCOEL 

N1=M+N-2 

TINV=TINVER(N1,ALPHA ) 

TINV2=-TINV 

JJ=1 

DO 220 NONCEN=MINDEL,MAXDEL,INCDEL 
TCDF=NONCDF(N1,NONCEN,TINV2 ) 

POW=TCOF 

NN( JJ )=NONCEN 

POWER( JJ )=POW 

JJ=JJ+1 

CONT INUE 

JJ=JJ-1 

PRINT, 'NN POWER' 

DO 230 I=1,JJ 

WRITE(6,24¢0) NN(I),POWER(T ) 

FORMAT('O', 'NN=‘,F8.354X, 'POWER=',F10.5) 
CONTINUE 

CALL PLOTT(NN,POWER >JJ>,0) 

PRINT*¥,'DO YOU WANT TO PLOT AGAIN (Y/N)?' 
READ (5,10) ANS 
IF(CANS.EQ.'Y').OR.(ANS.EQ.'Y')) GO TO 100 
GO TO 80 


TWO-SAMPLE TWO-SIOEO TEST, HA : MUX =/ MUY. 


PRINT*,'DO YOU WANT TO USE TWO-SIOED TEST(Y/N)?' 
READ(5,10) ANS 

IF((ANS.EQ.'N' ).OR.(ANS.EQ.'N')) GO TO 80 
PRINT¥,* MAXIMUM NONCENTRALITY VALUE 7? ' 
READ (55%) MAXDEL 

PRINT*,' MINIMUM NONCENTRALITY VALUE?‘ 
REAO (5,*) MINDEL 

PRINT*¥,*' INCREMENT NONCENTRALITY VALUE?‘ 
READ (5,x* )JINCOEL 

N1=M+N-2 

ALPHA1=ALPHA/2. 

TINV=TINVER(N1,ALPHAI ) 

JJ=1 

0O 250 NONCEN=MINOEL ,MAXOEL,INCOEL 
TCDF1=NONCDF(N1 »,NONCEN,TINV ) 
POW1=1.-TCDF1 

TINV1=-TINV 

TCDF2=NONCOF(N1,NONCEN,TINV1 ) 

POW2=TCDF2 

POW=POWN1+POW2 

NN( JJ J=NONCEN 

POWER( JJ J=POW 

JJ=JJ+1 

CONTINUE 

JJ=JJ-1 

PRINT, 'NN POWER‘ 

00 260 I=1,JJ 

WRITE(6,270) NN(I),POWER(T) 

FORMAT('O', 'NN=',F8.3,4X,'POWER=',F10.5) 
CONTINUE 

CALL PLOTT(NN,POWER ,JJ,0) 

PRINT*,'O0O YOU WANT TO PLOT AGAIN (Y/N)?' 
REAO (5,10) ANS 
IF((ANS.EQ.'Y').OR.CANS.EQ.'Y')) GO TO 100 
GO TO 80 

STOP 

ENO 


APPROXIMATION TO THE T- INVERSE OISTRIBUTION. 


CALCULATE T-INVERSE GIVEN OOF ALPHA. 
Q( T-ALPHA [| OOF )= ALPHA. 


FUNCTION TINVER(N,ALPHA ) 

REAL ALPHA,NUMER »,OENUM,LAM>XP,G1,G2,G3,X»TINVER 
INTEGER N 

DATA C0,C1,C2/2.515517,0.802853,0.010328/ 

OATA 01,02,03/1.432788,0.189269,0.001308/ 


CALCULATE NORMAL QUANTILE. 
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C1=(1.-€1.7(4¢.*NF ) ) XTSTAR-TLAM 


T=SQRT( ALOG(1./( ALPHA**2.))) 
NUME n=CO+C1*T+C2*T**2 

DENUM=1. +D1¥T+D2*T##2. +D3*TH*3. 
XP =T -(NUMER/DENUM ) 


CALCULATE THE INVERSE T-DISTRIBUTION. 


X=XP 


G1=.25%( X**5 +X ) 
G2=(1.796. J¥(5. XXXH54+16 . HXXHS+3 . €X ) 


G3=(1./7384¢. )*(3 . XXHK 7419, XXHHSZ417. XXX*5-15.*X ) 
G4=(1.792160. )¥0 79. XX*XH94+776 . KXHH74+1482 . XXH*X5-1592 . HXHHZ-9G5, €X ) 
TP=X+(G1/N )+(G2/(N¥N ) 40G3/0N%*3 ) 4+ 064/( N**4 ) ) 


TINVER=TP 


RETURN 


END 


THE OBJECTIVE OF THIS PROGRAM IS TO CALCULATE NON-CENTRAL T 


DISTRIBUTION (CDF ) 


THE CALCULATION METHOD IS NORMAL APPROXIMATION. 


HK KKK KKHKK HK KK KHKKEKHKHK VARIABLE DECLARATION  3¢3€3636363636363636 26 269 IEH HE HEHEHE HE EIEN 


FUNCTION NONCDF(NF »TLAM,TSTAR) 
INTEGER NF 
REAL TSTAR,TLAM,TX, NONCDF ,CDF 


DATA Al ,A2,A35P/.4361836,-.1201676> .9372980,.33267/ 
CONST 


= 1.0 7 SQRT (2.0 * 3.1415927) 


CALCULATE THE NORMAL APPROXIMATION . 


C2=SQRT(1.+( TSTAR¥*2/02.*NF ))) 
TX=C1/C2 

CDF = 0.0 

IF(TX .£Q.0) THEN 

CDF = 0.5 

ELSE IF(TX.GT.0) THEN 


T=1./01.+PTX ) 
XSQU=(TX*TX )*.5 
IF(XSQU.LT.17¢) THEN 

ZPDF =CONST XEXP( -XSQU ) 


CDF=1-ZPDF*( ALXT+AZXTHXZ+A3*T HS ) 

ELSE 

TX=-TX 

T=1./701.+4P*TX ) 

XSQU=( TX*¥TX)*.5 

IF(XSQU.LT.17¢) THEN 
ZPDF=CONST*EXP( -xXSQU } 


CDF=ZPDFHCUALET tAZHT HH Zt ASHT H3 ) 


END IF 


NONCDF=CDF 


RETURN 


END 
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